###############################################################################
############# Script for replication - APPENDIX DATAVERSE  ####################
###############################################################################
# R version 4.3.1
# (!) run main script(Script_APSR_main.R) before running the manuscript for reproducing the online appendix

######### APPENDIX D: Additional analyses #########

##### Check for country differences #####
## Reproduce: 
##      - Table D.1 (Appendix)  
##      - Figure D.1 (Appendix)
##      - Table D.2 (Appendix)
##      - Figure D.2 (Appendix)
##      - Table D.3 (Appendix)
##      - Figure D.3 (Appendix)
##      - Table D.4 (Appendix)
##
## Output: 
##      Table D.1: Share of respondents having chosen a scenario with certain attributes across countries 
##	    Figure D.1: Differences (AMCE) across countries, pairwise comparison
##	    Table D.2: Differences (AMCE) across countries, pairwise comparison
##	    Figure D.2: Differences (Marginal Means) across countries, pairwise comparison
##	    Table D.3: Differences (Marginal Means) across countries, pairwise comparison
##	    Figure D.3: AMCE country comparison and difference effects
##	    Table D. 4: AMCE country comparison and difference effects

# Differences across countries
cj_all %>% 
  group_by(issue, country) %>% 
  summarize(avg = mean(chosen))

cj_all %>% 
  group_by(initiative, country) %>% 
  summarize(avg = mean(chosen))

cj_all %>% 
  group_by(recruitment, country) %>% 
  summarize(avg = mean(chosen))

cj_all %>% 
  group_by(size, country) %>% 
  summarize(avg = mean(chosen))

cj_all %>% 
  group_by(composition, country) %>% 
  summarize(avg = mean(chosen))

cj_all %>% 
  group_by(majority, country) %>% 
  summarize(avg = mean(chosen))

cj_all %>% 
  group_by(aim, country) %>% 
  summarize(avg = mean(chosen))

cj_all %>% 
  group_by(authorization, country) %>% 
  summarize(avg = mean(chosen))

cj_all %>% 
  group_by(outcome.fav, country) %>% 
  summarize(avg = mean(chosen))

# make three data sets to compare effects pairwise 
FinIre <- cj_all |> filter(country != "USA")
FinUsa <- cj_all |> filter(country != "Ireland")
UsaIre <- cj_all |> filter(country != "Finland")

# Differences AMCE vs MM for Ireland - Finland 
prop.table(table(FinIre$country))
FinIre$country<- dplyr::recode(FinIre$country, Finland="Finland", Ireland="Ireland", USA="Ireland")
prop.table(table(FinIre$country))

diff_amce_FinIre_choice <- cj(FinIre, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)
diff_mm_FinIre_choice <- cj(FinIre, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names,feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)

# Differencs AMCE vs MM for USA - Finland 
prop.table(table(FinUsa$country))
FinUsa$country<- dplyr::recode(FinUsa$country, Finland="Finland", Ireland="USA", USA="USA")
prop.table(table(FinUsa$country))

diff_amce_FinUsa_choice <- cj(FinUsa, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)
diff_mm_FinUsa_choice <- cj(FinUsa, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names,feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)

# Differences AMCE vs MM for USA - Ireland 
prop.table(table(UsaIre$country))
UsaIre$country<- dplyr::recode(UsaIre$country, Finland="Ireland", Ireland="Ireland", USA="USA")
prop.table(table(UsaIre$country))

diff_amce_UsaIre_choice <- cj(UsaIre, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)
diff_mm_UsaIre_choice <- cj(UsaIre, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names,feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)

# Figure D. 1: Differences (AMCE) across countries, pairwise comparison
diff_amce_FinIre_choice$Estimate <- "Ireland-Finland"
diff_amce_FinUsa_choice$Estimate <- "USA-Finland"
diff_amce_UsaIre_choice$Estimate <- "USA-Ireland"

combined_amce <- rbind(diff_amce_FinIre_choice, diff_amce_FinUsa_choice, diff_amce_UsaIre_choice)

plot(combined_amce, feature_headers = FALSE) + 
  ggplot2::facet_wrap(~Estimate, ncol = 3) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 12),
        axis.title = element_blank(),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_blank(),
        strip.text.x = element_text(colour = "#091926", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = c("#515151","#515151","#515151"))

# Table D. 2: Differences (AMCE) across countries, pairwise comparison
C1<- head(combined_amce[c("level", "estimate", "std.error")], 30L)
knitr::kable(C1, format = "markdown")

# Figure D. 2: Differences (Marginal Means) across countries, pairwise comparison 
diff_mm_FinIre_choice$Estimate <- "Ireland-Finland"
diff_mm_FinUsa_choice$Estimate <- "USA-Finland"
diff_mm_UsaIre_choice$Estimate <- "USA-Ireland"

combined_mm <- rbind(diff_mm_FinIre_choice, diff_mm_FinUsa_choice, diff_mm_UsaIre_choice)

plot(combined_mm, feature_headers = FALSE) + 
  ggplot2::facet_wrap(~Estimate, ncol = 3) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 12),
        axis.title = element_blank(),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_blank(),
        strip.text.x = element_text(colour = "#091926", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = c("#515151","#515151","#515151","#515151","#515151","#515151","#515151","#515151","#515151"))

# Table D. 3: Differences (Marginal Means) across countries, pairwise comparison
C2<- head(combined_mm[c("level", "estimate", "std.error")], 57L)
knitr::kable(C2, format = "markdown")

# Figure D. 3: AMCE country comparison and difference effects
amce_by_FinIre_choice <- cj(FinIre, Mchoice, id = ~ID,estimate = "amce",feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)
diff_amce_FinIre_choice <- cj(FinIre, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)

amce_FinIre_choice_diff <- rbind(amce_by_FinIre_choice, diff_amce_FinIre_choice)

plot(amce_FinIre_choice_diff, feature_headers = FALSE) +
  ggplot2::facet_grid(feature~BY, scales = "free_y", space="free", switch = "y") +
  theme(strip.background.x = element_blank(),
        strip.text.x = element_text(colour = "#091926", size = 14),
        strip.text.y = element_blank(),
        panel.spacing.y = unit(1, "lines"),
        panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#091926", size = 12),
        axis.title = element_text(colour = "#091926", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = c("#515151", "#515151"))

amce_by_FinUsa_choice <- cj(FinUsa, Mchoice, id = ~ID,estimate = "amce",feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)
diff_amce_FinUsa_choice <- cj(FinUsa, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)

amce_FinUsa_choice_diff <- rbind(amce_by_FinUsa_choice, diff_amce_FinUsa_choice)

plot(amce_FinUsa_choice_diff, feature_headers = FALSE) +
  ggplot2::facet_grid(feature~BY, scales = "free_y", space="free", switch = "y") +
  theme(strip.background.x = element_blank(),
        strip.text.x = element_text(colour = "#091926", size = 14),
        strip.text.y = element_blank(),
        panel.spacing.y = unit(1, "lines"),
        panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#091926", size = 12),
        axis.title = element_text(colour = "#091926", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = c("#515151", "#515151"))

amce_by_UsaIre_choice <- cj(UsaIre, Mchoice, id = ~ID,estimate = "amce",feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)
diff_amce_UsaIre_choice <- cj(UsaIre, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~country)

amce_UsaIre_choice_diff <- rbind(amce_by_UsaIre_choice, diff_amce_UsaIre_choice)

plot(amce_UsaIre_choice_diff, feature_headers = FALSE) +
  ggplot2::facet_grid(feature~BY, scales = "free_y", space="free", switch = "y") +
  theme(strip.background.x = element_blank(),
        strip.text.x = element_text(colour = "#091926", size = 14),
        strip.text.y = element_blank(),
        panel.spacing.y = unit(1, "lines"),
        panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#091926", size = 12),
        axis.title = element_text(colour = "#091926", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = c("#515151", "#515151"))

# Table D. 4: AMCE country comparison and difference effects
C3a <- head(amce_FinIre_choice_diff[c("level", "estimate", "std.error")], 48L)
knitr::kable(C3a, format = "markdown")

C3b <- head(amce_FinUsa_choice_diff [c("level", "estimate", "std.error")], 48L)
knitr::kable(C3b, format = "markdown")

C3c <- head(amce_UsaIre_choice_diff [c("level", "estimate", "std.error")], 48L)
knitr::kable(C3c, format = "markdown")

##### A) All countries  #####
#### Familiarity and experience ####

## Reproduce: 
##      - Figure D.4 (Appendix)
##      - Table D.5 (Appendix)
##      - Figure D.5 (Appendix)
##      - Table D.6 (Appendix)
##
## Output: 
##      Figure D.4: Marginal means for FAMILIARITY, PARTICIPATED and EXPERIENCES 
##	    Table D.5: Marginal means for FAMILIARITY, PARTICIPATED and EXPERIENCES
##	    Figure D.5: AMCE for respondents who have participated (1) and who have not participated (0)
##	    Table D.6: AMCE for respondents who have participated (1) and who have not participated (0)
##	    

# Figure D. 4: Marginal means for FAMILIARITY, PARTICIPATED and EXPERIENCES
diff_mm_exp_choice <- cj(cj_experience, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~familiar)
diff_mm_participated_choice <- cj(cj_experience, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~participated)
diff_mm_expdirect_choice <- cj(cj_expdirect, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~exsplit)

diff_experience_estimates_choice_mm <-bind_rows(
  diff_mm_exp_choice %>% mutate(model = "Familiarity"),
  diff_mm_participated_choice %>% mutate(model = "Participated"),
  diff_mm_expdirect_choice %>% mutate(model = "Positive experience"))

diff_experience_choice_mm <- diff_experience_estimates_choice_mm %>% select(BY:level, estimate:model)
diff_experience_choice_mm$model <- factor(diff_experience_choice_mm$model, levels = c("Familiarity", "Participated", "Positive experience"))

plot(diff_experience_choice_mm, feature_headers = FALSE, xlab = "Estimated Difference (Marginal Means)") + 
  ggplot2::facet_wrap(~model, ncol = 3) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = rep("#515151", 9))

# Table D.5: Marginal means for FAMILIARITY, PARTICIPATED and EXPERIENCES
head(diff_experience_choice_mm[c("level", "estimate", "std.error")], 57L)

# Figure D.5: AMCE for respondents who have participated (1) and who have not participated (0)
# estimate two conjoint models for "having participated" vs. "rest"
amce_by_participated <- cj(cj_experience, Mchoice, id = ~ID, estimate = "amce", feature_labels = feature_names, feature_order = c("outcome.fav",  "authorization", "aim", "majority", "composition",  "size", "recruitment", "initiative", "issue"), by = ~participated)

plot(amce_by_participated, group = "participated", feature_headers = FALSE)+
  ggplot2::facet_wrap(~BY, ncol = 2L) +
  theme(panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 12),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none",
        panel.spacing.y = unit(1, "lines"),
        strip.background= element_rect(color="white", size = 2, fill = "gray97"), 
        strip.text.y.left = element_text(angle = 0),
        strip.placement = "outside",
        strip.text.y = element_text(colour = "#515151", size = 10, face = "bold")) +
  theme(plot.title = element_text(colour =  "#515151", size = 14)) +
  geom_point(size=2) +
  scale_color_manual(values = c("#515151", "#515151", "#515151"))

# Table D.6:AMCE for respondents who have participated (1) and who have not participated (0)
C5 <- head(amce_by_participated [c("level", "estimate", "std.error")], 38L)
knitr::kable(C5, format = "markdown")

#### Trust ####
## Reproduce: 
##      - Figure D.6 (Appendix)
##      - Table D.7 (Appendix)
##      - Figure D.7 (Appendix)
##      - Table D.8 (Appendix)
##
## Output: 
##      Figure D.6: Marginal means LOW POLITICAL TRUST (top row) and HIGH TRUST IN CITIZENS 
##	    Table D.7: Marginal means LOW POLITICAL TRUST (top row) and HIGH TRUST IN CITIZENS
##	    Figure D.7: Marginal means for respondents with LOW POLITICAL trust but HIGH TRUST IN CITIZENS
##	    Table D.8: Marginal means for respondents with LOW POLITICAL trust but HIGH TRUST IN CITIZENS
##
# Figure D.6: Marginal means LOW POLITICAL TRUST (top row) and HIGH TRUST IN CITIZENS 
diff_mm_politics <- cj(cj_trust_politics, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
diff_mm_politicians <- cj(cj_trust_politicians, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
diff_mm_dmp <- cj(cj_trust, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~dmpsplit)
diff_mm_citizens <- cj(cj_trust, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~citizensplit)

diff_disent_estimates_trust_mm <-bind_rows(
  diff_mm_politics %>% mutate(model = "Political institutions"),
  diff_mm_politicians %>% mutate(model = "Politicians"),
  diff_mm_dmp %>% mutate(model = "Minipublics"),
  diff_mm_citizens %>% mutate(model = "Citizens"))

diff_disent_trust_mm <- diff_disent_estimates_trust_mm %>% select(BY:level, estimate:model)
diff_disent_trust_mm$model <- factor(diff_disent_trust_mm$model, levels = c("Political institutions", "Politicians", "Minipublics", "Citizens"))

plot(diff_disent_trust_mm, feature_headers = FALSE, xlab = "Estimated Difference (Marginal Means)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = rep("#515151", 9))

# Table D.7: Marginal means LOW POLITICAL TRUST (top row) and HIGH TRUST IN CITIZENS
head(diff_disent_trust_mm[c("level", "estimate", "std.error")], 76L)

# Figure D.7: Marginal means for respondents with LOW POLITICAL trust but HIGH TRUST IN CITIZENS
diff_mm_select_all <- cj(cj_trust_NA, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~select_trust)
diff_mm_select_fin <- cj(cj_fin_NA, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~select_trust)
diff_mm_select_ire<- cj(cj_ire_NA, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~select_trust)
diff_mm_select_usa <- cj(cj_usa_NA, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~select_trust)

diff_disent_estimates_select_mm <-bind_rows(
  diff_mm_select_all %>% mutate(model = "All countries"),
  diff_mm_select_fin %>% mutate(model = "Finland"),
  diff_mm_select_ire %>% mutate(model = "Ireland"),
  diff_mm_select_usa %>% mutate(model = "USA"))

diff_disent_select_mm <- diff_disent_estimates_select_mm %>% select(BY:level, estimate:model)
diff_disent_select_mm$model <- factor(diff_disent_select_mm$model, levels = c("All countries", "Finland", "Ireland", "USA"))

plot(diff_disent_select_mm, feature_headers = FALSE, xlab = "Estimated Difference (Marginal Means)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = rep("#515151", 9))

# Table D.8: Marginal means for respondents with LOW POLITICAL trust but HIGH TRUST IN CITIZENS
head(diff_disent_select_mm[c("level", "estimate", "std.error")], 76L)

##### B) Finland #####
#### Familiarity and experience ####
## Reproduce: 
##      - Figure D.8 (Appendix)
##      - Table D.9 (Appendix)
##      - Figure D.9 (Appendix)
##      - Table D.10 (Appendix)
##
## Output: 
##      Figure D.8: AMCE of minipublic attributes on support contingent on familiarity and experiences - Finnish sample
##	    Table D.9: AMCE of minipublic attributes on support contingent on familiarity and experiences - Finnish sample
##	    Figure D.9: MARGINAL MEANS for familiarity and experiences – Finnish sample
##	    Table D.10: MARGINAL MEANS for familiarity and experiences – Finnish sample
##
# Dummy (familiar vs. unfamiliar): at least know well vs. unfamiliar and heard about
experiences_fin$familiar <- ifelse(experiences_fin$experiences > 2, 1, 0)
experiences_fin$familiar <- as.factor(experiences_fin$familiar)
summary(experiences_fin$familiar)

# Dummy (participated vs. not)
experiences_fin$participated <- ifelse(experiences_fin$experiences == 4, 1, 0)
experiences_fin$participated <- as.factor(experiences_fin$participated)

summary(experiences_fin$participated)

# Direction of experience
describe(experiences_fin$experiences_direct)
fin_expdirect<- experiences_fin %>% select(respondent_id:chosen, outcome.fav, experiences, experiences_direct, retro, ID)
fin_expdirect <-  fin_expdirect %>% filter(experiences != 1)
fin_expdirect$exsplit<-dicho(fin_expdirect$experiences_direct, dich.by = "median", as.num = TRUE)
fin_expdirect$exsplit <- factor(fin_expdirect$exsplit, labels = c("negative experiences", "positive experiences"))

summary(fin_expdirect$exsplit)

cj_anova(fin_expdirect, Mchoice, by = ~exsplit)

t.test(experiences_fin$retro ~ experiences_fin$familiar, var.equal = TRUE, alternative = "two.sided") # significant differences (familiar: mean=4.83; unfamiliar: mean=4.37)
cohensD(experiences_fin$retro ~ experiences_fin$familiar) # Cohens d = 0.36 (small effect)

t.test(fin_expdirect$retro ~ fin_expdirect$exsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (positive: mean=5.16; negative: mean=4.22)
cohensD(fin_expdirect$retro ~ fin_expdirect$exsplit) # Cohens d = 0.87 (large effect)

# Figure D. 8: AMCE of minipublic attributes on support contingent on familiarity and experiences - Finnish sample
diff_amce_exp_choice_fin <- cj(experiences_fin, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~familiar)
diff_amce_participated_choice_fin <- cj(experiences_fin, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~participated)
diff_amce_expdirect_choice_fin <- cj(fin_expdirect, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~exsplit)

diff_experience_estimates_choice_fin <-bind_rows(
  diff_amce_exp_choice_fin %>% mutate(model = "Familiarity"),
  diff_amce_participated_choice_fin %>% mutate(model = "Participated"),
  diff_amce_expdirect_choice_fin %>% mutate(model = "Positive experience"))

diff_experience_choice_fin <- diff_experience_estimates_choice_fin %>% select(BY:level, estimate:model)
diff_experience_choice_fin$model <- factor(diff_experience_choice_fin$model, levels = c("Familiarity", "Participated", "Positive experience"))

plot(diff_experience_choice_fin, feature_headers = FALSE, xlab = "Estimated Difference (AMCE)") + 
  ggplot2::facet_wrap(~model, ncol = 3) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("Finland") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 9: AMCE of minipublic attributes on support contingent on familiarity and experiences - Finnish sample
head(diff_experience_choice_fin[c("level", "estimate", "std.error")], 30L)

# Figure D. 9: MARGINAL MEANS for familiarity and experiences – Finnish sample
diff_mm_exp_choice_fin <- cj(experiences_fin, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~familiar)
diff_mm_participated_choice_fin <- cj(experiences_fin, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~participated)
diff_mm_expdirect_choice_fin <- cj(fin_expdirect, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~exsplit)

diff_experience_estimates_choice_fin_mm <-bind_rows(
  diff_mm_exp_choice_fin %>% mutate(model = "Familiarity"),
  diff_mm_participated_choice_fin %>% mutate(model = "Participated"),
  diff_mm_expdirect_choice_fin %>% mutate(model = "Positive experience"))

diff_experience_choice_fin_mm <- diff_experience_estimates_choice_fin_mm %>% select(BY:level, estimate:model)
diff_experience_choice_fin_mm$model <- factor(diff_experience_choice_fin_mm$model, levels = c("Familiarity", "Participated", "Positive experience"))

plot(diff_experience_choice_fin_mm, feature_headers = FALSE, xlab = "Estimated Difference (Marginal Means)") + 
  ggplot2::facet_wrap(~model, ncol = 3) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("Finland") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 10: MARGINAL MEANS for familiarity and experiences – Finnish sample
head(diff_experience_choice_fin_mm[c("level", "estimate", "std.error")], 57L)

#### Trust ####
## Reproduce: 
##      - Figure D.10 (Appendix)
##      - Table D.11 (Appendix)
##      - Figure D.11 (Appendix)
##      - Table D.12 (Appendix)
##
## Output: 
##      Figure D.10: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citi-zens (bottom row) – Finnish sample
##	    Table D.11: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citizens (bottom row) – Finnish sample
##	    Figure D.11: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – Finnish sample
##	    Table D.12: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – Finnish sample
##
# Trust in political institutions: index parliament and government
alpha(subset(cj_fin, select = c(trust_parliament, trust_government)), check.keys = TRUE) # alpha .9
cj_fin$trust_politics <- rowMeans(subset(cj_fin, select = c(trust_parliament, trust_government)), na.rm = TRUE)
describe(cj_fin$trust_politics) # mean: 4.17, median: 4.5
cj_fin_politics <- cj_fin %>% filter(trust_politics != "NaN")
cj_fin_politics$trustsplit<-dicho(cj_fin_politics$trust_politics, dich.by = "median", as.num = TRUE)
class(cj_fin_politics$trustsplit)
cj_fin_politics$trustsplit <- factor(cj_fin_politics$trustsplit, labels = c("Low trust", "High trust"))
cj_fin_politics$trustsplit <- relevel(cj_fin_politics$trustsplit, ref = "High trust")
cj_anova(cj_fin_politics, Mchoice, by = ~trustsplit)

summary(cj_fin_politics$trustsplit)

# Trust in politicians
describe(cj_fin$trust_politicians) # mean: 3.71, median: 4
fin_politicians_NA <- is.na(select(cj_fin, trust_politicians))
cj_fin_politicians <- filter(cj_fin, fin_politicians_NA == 0)
describe(cj_fin_politicians$trust_politicians) #mean 3.71 median 4
cj_fin_politicians$trustsplit<-dicho(cj_fin_politicians$trust_politicians, dich.by = "median", as.num = TRUE)
class(cj_fin_politicians$trustsplit)
cj_fin_politicians$trustsplit <- factor(cj_fin_politicians$trustsplit, labels = c("Low trust", "High trust"))
cj_fin_politicians$trustsplit <- relevel(cj_fin_politicians$trustsplit, ref = "High trust")
cj_anova(cj_fin_politicians, Mchoice, by = ~trustsplit)

summary(cj_fin_politicians$trustsplit)

# Trust in minipublics
describe(cj_fin$trust_dmp) # mean: 4.2, median: 4
cj_fin$dmpsplit<-dicho(cj_fin$trust_dmp, dich.by = "median", as.num = TRUE)
class(cj_fin$dmpsplit)
cj_fin$dmpsplit <- factor(cj_fin$dmpsplit, labels = c("Low trust", "High trust"))
cj_anova(cj_fin, Mchoice, by = ~dmpsplit)

summary(cj_fin$dmpsplit)

# Trust in other citizens
alpha(subset(cj_fin, select = c(trust_citizens1, trust_citizens2)), check.keys = TRUE) # alpha .86
cj_fin$trust_citizens <- rowMeans(subset(cj_fin, select = c(trust_citizens1, trust_citizens2)), na.rm = TRUE)
describe(cj_fin$trust_citizens) #mean 4.3 median 4.5
cj_fin$citizensplit<-dicho(cj_fin$trust_citizens, dich.by = "median", as.num = TRUE)
class(cj_fin$citizensplit)
cj_fin$citizensplit <- factor(cj_fin$citizensplit, labels = c("Low trust", "High trust"))
cj_anova(cj_fin, Mchoice, by = ~citizensplit)

summary(cj_fin$citizensplit)

describeBy(cj_fin_politics$retro, cj_fin_politics$trustsplit)
describeBy(cj_fin_politicians$retro, cj_fin_politicians$trustsplit)
describeBy(cj_fin$retro, cj_fin$dmpsplit)
describeBy(cj_fin$retro, cj_fin$citizensplit)

t.test(cj_fin_politics$retro ~ cj_fin_politics$trustsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=4.73; low trust: mean=4.29)
cohensD(cj_fin_politics$retro ~ cj_fin_politics$trustsplit) # Cohens d = 0.35 (small effect)

t.test(cj_fin_politicians$retro ~ cj_fin_politicians$trustsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=4.74; low trust: mean=4.35)
cohensD(cj_fin_politicians$retro ~ cj_fin_politicians$trustsplit) # Cohens d = 0.29 (small effect)

t.test(cj_fin$retro ~ cj_fin$dmpsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.27; low trust: mean=3.93)
cohensD(cj_fin$retro ~ cj_fin$dmpsplit) # Cohens d = 1.21 (large effect)

t.test(cj_fin$retro ~ cj_fin$citizensplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.17; low trust: mean=4.03)
cohensD(cj_fin$retro ~ cj_fin$citizensplit) # Cohens d = 0.98 (large effect)

# Figure D. 10: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citi-zens (bottom row) – Finnish sample
fin_diff_amce_politics <- cj(cj_fin_politics, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
fin_diff_amce_politicians <- cj(cj_fin_politicians, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
fin_diff_amce_dmp <- cj(cj_fin, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~dmpsplit)
fin_diff_amce_citizens <- cj(cj_fin, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~citizensplit)

fin_diff_disent_estimates_trust <-bind_rows(
  fin_diff_amce_politics %>% mutate(model = "Political institutions"),
  fin_diff_amce_politicians %>% mutate(model = "Politicians"),
  fin_diff_amce_dmp %>% mutate(model = "Minipublics"),
  fin_diff_amce_citizens %>% mutate(model = "Citizens"))

fin_diff_disent_trust <- fin_diff_disent_estimates_trust %>% select(BY:level, estimate:model)
fin_diff_disent_trust$model <- factor(fin_diff_disent_trust$model, levels = c("Political institutions", "Politicians", "Minipublics", "Citizens"))

plot(fin_diff_disent_trust, feature_headers = FALSE, xlab = "Estimated Difference (AMCE)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("Finland") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 11: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citizens (bottom row) – Finnish sample
head(fin_diff_disent_trust[c("level", "estimate", "std.error")], 40L)

# Figure D. 11: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – Finnish 
fin_diff_mm_politics <- cj(cj_fin_politics, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
fin_diff_mm_politicians <- cj(cj_fin_politicians, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
fin_diff_mm_dmp <- cj(cj_fin, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~dmpsplit)
fin_diff_mm_citizens <- cj(cj_fin, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~citizensplit)

fin_diff_disent_estimates_trust_mm <-bind_rows(
  fin_diff_mm_politics %>% mutate(model = "Political institutions"),
  fin_diff_mm_politicians %>% mutate(model = "Politicians"),
  fin_diff_mm_dmp %>% mutate(model = "Minipublics"),
  fin_diff_mm_citizens %>% mutate(model = "Citizens"))

fin_diff_disent_trust_mm <- fin_diff_disent_estimates_trust_mm %>% select(BY:level, estimate:model)
fin_diff_disent_trust_mm$model <- factor(fin_diff_disent_trust_mm$model, levels = c("Political institutions", "Politicians", "Minipublics", "Citizens"))

plot(fin_diff_disent_trust_mm, feature_headers = FALSE, xlab = "Estimated Difference (Marginal Means)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("Finland") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 12: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – Finnish sample
head(fin_diff_disent_trust_mm[c("level", "estimate", "std.error")], 76L)

##### C) Ireland #####
#### Familiarity and experience ####
## Reproduce: 
##      - Figure D.12 (Appendix)
##      - Table D.13 (Appendix)
##      - Figure D.13 (Appendix)
##      - Table D.14 (Appendix)
##
## Output: 
##      Figure D.12: AMCE of minipublic attributes on support contingent on familiarity and experiences - Irish sample
##	    Table D.13: AMCE of minipublic attributes on support contingent on familiarity and experiences - Irish sample
##	    Figure D.13: MARGINAL MEANS for familiarity and experiences – Irish sample
##	    Table D.14: MARGINAL MEANS for familiarity and experiences – Irish sample
##
# Dummy (familiar vs. unfamiliar): at least know well vs. unfamiliar and heard about
experiences_ire$familiar <- ifelse(experiences_ire$experiences > 2, 1, 0)
experiences_ire$familiar <- as.factor(experiences_ire$familiar)

summary(experiences_ire$familiar)

# Dummy (participated vs. not)
experiences_ire$participated <- ifelse(experiences_ire$experiences == 4, 1, 0)
experiences_ire$participated <- as.factor(experiences_ire$participated)

summary(experiences_ire$participated)

# Direction of experiences (positive vs negative) 
describe(experiences_ire$experiences_direct)
ire_expdirect<- experiences_ire %>% select(respondent_id:chosen, outcome.fav, experiences, experiences_direct, retro, ID) 
ire_expdirect <- ire_expdirect %>% filter(experiences != 1)
ire_expdirect$exsplit<-dicho(ire_expdirect$experiences_direct, dich.by = "median", as.num = TRUE)
ire_expdirect$exsplit <- factor(ire_expdirect$exsplit, labels = c("negative experiences", "positive experiences"))

summary(ire_expdirect$exsplit)

cj_anova(ire_expdirect, Mchoice, by = ~exsplit)

t.test(experiences_ire$retro ~ experiences_ire$familiar, var.equal = TRUE, alternative = "two.sided") # significant differences (familiar: mean=5.31; unfamiliar: mean=4.97)
cohensD(experiences_ire$retro ~ experiences_ire$familiar) # Cohens d = 0.26 (small effect)

t.test(ire_expdirect$retro ~ ire_expdirect$exsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (positive: mean=6.06; negative: mean=4.83)
cohensD(ire_expdirect$retro ~ ire_expdirect$exsplit) # Cohens d = 1.02 (large effect)

# Figure D. 12: AMCE of minipublic attributes on support contingent on familiarity and experiences - Irish sample
diff_amce_exp_choice_ire <- cj(experiences_ire, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~familiar)
diff_amce_participated_choice_ire <- cj(experiences_ire, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~participated)
diff_amce_expdirect_choice_ire <- cj(ire_expdirect, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~exsplit)

diff_experience_estimates_choice_ire <-bind_rows(
  diff_amce_exp_choice_ire %>% mutate(model = "Familiarity"),
  diff_amce_participated_choice_ire %>% mutate(model = "Participated"),
  diff_amce_expdirect_choice_ire %>% mutate(model = "Positive experience"))

diff_experience_choice_ire <- diff_experience_estimates_choice_ire %>% select(BY:level, estimate:model)
diff_experience_choice_ire$model <- factor(diff_experience_choice_ire$model, levels = c("Familiarity", "Participated", "Positive experience"))

plot(diff_experience_choice_ire, feature_headers = FALSE, xlab = "Estimated Difference (AMCE)") + 
  ggplot2::facet_wrap(~model, ncol = 3) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("Ireland") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 13: AMCE of minipublic attributes on support contingent on familiarity and experiences - Irish sample
head(diff_experience_choice_ire[c("level", "estimate", "std.error")], 30L)

# Figure D. 13: MARGINAL MEANS for familiarity and experiences – Irish sample
diff_mm_exp_choice_ire <- cj(experiences_ire, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~familiar)
diff_mm_participated_choice_ire <- cj(experiences_ire, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~participated)
diff_mm_expdirect_choice_ire <- cj(ire_expdirect, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~exsplit)

diff_experience_estimates_choice_ire_mm <-bind_rows(
  diff_mm_exp_choice_ire %>% mutate(model = "Familiarity"),
  diff_mm_participated_choice_ire %>% mutate(model = "Participated"),
  diff_mm_expdirect_choice_ire %>% mutate(model = "Positive experience"))

diff_experience_choice_ire_mm <- diff_experience_estimates_choice_ire_mm %>% select(BY:level, estimate:model)
diff_experience_choice_ire_mm$model <- factor(diff_experience_choice_ire_mm$model, levels = c("Familiarity", "Participated", "Positive experience"))

plot(diff_experience_choice_ire_mm, feature_headers = FALSE, xlab = "Estimated Difference (Marginal Means)") + 
  ggplot2::facet_wrap(~model, ncol = 3) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("Ireland") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 14: MARGINAL MEANS for familiarity and experiences – Irish sample
head(diff_experience_choice_ire_mm[c("level", "estimate", "std.error")], 57L)

#### Trust ####
## Reproduce: 
##      - Figure D.14 (Appendix)
##      - Table D.15 (Appendix)
##      - Figure D.15 (Appendix)
##      - Table D.16 (Appendix)
##
## Output: 
##      Figure D.14: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citi-zens (bottom row) – Irish sample
##	    Table D.15: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citizens (bottom row) – Irish sample
##	    Figure D.15:  MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – Irish sample
##	    Table D.16: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – Irish sample
##
# Trust in political institutions: index parliament and government
alpha(subset(cj_ire, select = c(trust_parliament, trust_government)), check.keys = TRUE) # alpha .88
cj_ire$trust_politics <- rowMeans(subset(cj_ire, select = c(trust_parliament, trust_government)), na.rm = TRUE)
describe(cj_ire$trust_politics) # mean: 4 median: 4
cj_ire_politics <- cj_ire %>% filter(trust_politics != "NaN")
cj_ire_politics$trustsplit<-dicho(cj_ire_politics$trust_politics, dich.by = "median", as.num = TRUE)
class(cj_ire_politics$trustsplit)
cj_ire_politics$trustsplit <- factor(cj_ire_politics$trustsplit, labels = c("Low trust", "High trust"))
cj_ire_politics$trustsplit <- relevel(cj_ire_politics$trustsplit, ref = "High trust")
cj_anova(cj_ire_politics, Mchoice, by = ~trustsplit)

summary(cj_ire_politics$trustsplit)

# Trust in politicians
ire_politicians_NA <- is.na(select(cj_ire, trust_politicians))
cj_ire_politicians <- filter(cj_ire, ire_politicians_NA == 0)
describe(cj_ire_politicians$trust_politicians) #mean 3.45 median 4
cj_ire_politicians$trustsplit<-dicho(cj_ire_politicians$trust_politicians, dich.by = "median", as.num = TRUE)
class(cj_ire_politicians$trustsplit)
cj_ire_politicians$trustsplit <- factor(cj_ire_politicians$trustsplit, labels = c("Low trust", "High trust"))
cj_ire_politicians$trustsplit <- relevel(cj_ire_politicians$trustsplit, ref = "High trust")

summary(cj_ire_politicians$trustsplit)

cj_anova(cj_ire_politicians, Mchoice, by = ~trustsplit)

# Trust in minipublics
describe(cj_ire$trust_dmp) # mean: 4.81, mediam: 5
cj_ire$dmpsplit<-dicho(cj_ire$trust_dmp, dich.by = "median", as.num = TRUE)
class(cj_ire$dmpsplit)
cj_ire$dmpsplit <- factor(cj_ire$dmpsplit, labels = c("Low trust", "High trust"))

summary(cj_ire$dmpsplit)

cj_anova(cj_ire, Mchoice, by = ~dmpsplit)

# Trust in other citizens
alpha(subset(cj_ire, select = c(trust_citizens1, trust_citizens2)), check.keys = TRUE) # alpha .84
cj_ire$trust_citizens <- rowMeans(subset(cj_ire, select = c(trust_citizens1, trust_citizens2)), na.rm = TRUE)
describe(cj_ire$trust_citizens) # mean: 4.61, median: 4.5
cj_ire$citizensplit<-dicho(cj_ire$trust_citizens, dich.by = "median", as.num = TRUE)
class(cj_ire$citizensplit)
cj_ire$citizensplit <- factor(cj_ire$citizensplit, labels = c("Low trust", "High trust"))

summary(cj_ire$citizensplit)

cj_anova(cj_ire, Mchoice, by = ~citizensplit)

describeBy(cj_ire_politics$retro, cj_ire_politics$trustsplit)
describeBy(cj_ire_politicians$retro, cj_ire_politicians$trustsplit)
describeBy(cj_ire$retro, cj_ire$dmpsplit)
describeBy(cj_ire$retro, cj_ire$citizensplit)

t.test(cj_ire_politics$retro ~ cj_ire_politics$trustsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.26; low trust: mean=4.87)
cohensD(cj_ire_politics$retro ~ cj_ire_politics$trustsplit) # Cohens d = 0.31 (small effect)

t.test(cj_ire_politicians$retro ~ cj_ire_politicians$trustsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.27; low trust: mean=4.97)
cohensD(cj_ire_politicians$retro ~ cj_ire_politicians$trustsplit) # Cohens d = 0.23 (small effect)

t.test(cj_ire$retro ~ cj_ire$dmpsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=6.08; low trust: mean=4.64)
cohensD(cj_ire$retro ~ cj_ire$dmpsplit) # Cohens d = 1.26 (large effect)

t.test(cj_ire$retro ~ cj_ire$citizensplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.65; low trust: mean=4.48)
cohensD(cj_ire$retro ~ cj_ire$citizensplit) # Cohens d = 1.01 (large effect)

# Figure D. 14: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citi-zens (bottom row) – Irish sample
ire_diff_amce_politics <- cj(cj_ire_politics, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
ire_diff_amce_politicians <- cj(cj_ire_politicians, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
ire_diff_amce_dmp <- cj(cj_ire, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~dmpsplit)
ire_diff_amce_citizens <- cj(cj_ire, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~citizensplit)

ire_diff_disent_estimates_trust <-bind_rows(
  ire_diff_amce_politics %>% mutate(model = "Political institutions"),
  ire_diff_amce_politicians %>% mutate(model = "Politicians"),
  ire_diff_amce_dmp %>% mutate(model = "Minipublics"),
  ire_diff_amce_citizens %>% mutate(model = "Citizens"))

ire_diff_disent_trust <- ire_diff_disent_estimates_trust %>% select(BY:level, estimate:model)
ire_diff_disent_trust$model <- factor(ire_diff_disent_trust$model, levels = c("Political institutions", "Politicians", "Minipublics", "Citizens"))

plot(ire_diff_disent_trust, feature_headers = FALSE, xlab = "Estimated Difference (AMCE)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("Ireland") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 15: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citizens (bottom row) – Irish sample
head(ire_diff_disent_trust[c("level", "estimate", "std.error")], 40L)

# Figure D. 15:  MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – Irish sample
ire_diff_mm_politics <- cj(cj_ire_politics, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
ire_diff_mm_politicians <- cj(cj_ire_politicians, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
ire_diff_mm_dmp <- cj(cj_ire, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~dmpsplit)
ire_diff_mm_citizens <- cj(cj_ire, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~citizensplit)

ire_diff_disent_estimates_trust_mm <-bind_rows(
  ire_diff_mm_politics %>% mutate(model = "Political institutions"),
  ire_diff_mm_politicians %>% mutate(model = "Politicians"),
  ire_diff_mm_dmp %>% mutate(model = "Minipublics"),
  ire_diff_mm_citizens %>% mutate(model = "Citizens"))

ire_diff_disent_trust_mm <- ire_diff_disent_estimates_trust_mm %>% select(BY:level, estimate:model)
ire_diff_disent_trust_mm$model <- factor(ire_diff_disent_trust_mm$model, levels = c("Political institutions", "Politicians", "Minipublics", "Citizens"))

plot(ire_diff_disent_trust_mm, feature_headers = FALSE, xlab = "Estimated Difference (Marginal Means)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("Ireland") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 16: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – Irish sample
head(ire_diff_disent_trust_mm[c("level", "estimate", "std.error")], 76L)

##### D) USA #####
#### Familiarity and experience ####
## Reproduce: 
##      - Figure D.16 (Appendix)
##      - Table D.17 (Appendix)
##      - Figure D.17 (Appendix)
##      - Table D.18 (Appendix)
##
## Output: 
##      Figure D.16: AMCE of minipublic attributes on support contingent on familiarity and experiences - US sample
##	    Table D.17: AMCE of minipublic attributes on support contingent on familiarity and experiences - US sample
##	    Figure D.17: MARGINAL MEANS for familiarity and experiences – US sample
##	    Table D.18: MARGINAL MEANS for familiarity and experiences – US sample
##
# Dummy (familiar vs. unfamiliar): at least know well vs. unfamiliar and heard about
experiences_usa$familiar <- ifelse(experiences_usa$experiences > 2, 1, 0)
experiences_usa$familiar <- as.factor(experiences_usa$familiar)

summary(experiences_usa$familiar)

# Dummy (participated vs. not)
experiences_usa$participated <- ifelse(experiences_usa$experiences == 4, 1, 0)
experiences_usa$participated <- as.factor(experiences_usa$participated)

summary(experiences_usa$participated)

# Direction of experiences (positive vs negative)
describe(experiences_usa$experiences_direct)
usa_expdirect<- experiences_usa %>% select(respondent_id:chosen, outcome.fav, experiences, experiences_direct, retro, ID)
usa_expdirect <- usa_expdirect %>% filter(experiences != 1)
usa_expdirect$exsplit<-dicho(usa_expdirect$experiences_direct, dich.by = "median", as.num = TRUE)
usa_expdirect$exsplit <- factor(usa_expdirect$exsplit, labels = c("negative experiences", "positive experiences"))

summary(usa_expdirect$exsplit)

cj_anova(usa_expdirect, Mchoice, by = ~exsplit)

t.test(experiences_usa$retro ~ experiences_usa$familiar, var.equal = TRUE, alternative = "two.sided") # significant differences (familiar: mean=5.47; unfamiliar: mean=4.94)
cohensD(experiences_usa$retro ~ experiences_usa$familiar) # Cohens d = 0.39 (small effect)

t.test(usa_expdirect$retro ~ usa_expdirect$exsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (positive: mean=6.13; negative: mean=4.89)
cohensD(usa_expdirect$retro ~ usa_expdirect$exsplit) # Cohens d = 1.06 (large effect)

# Figure D. 16: AMCE of minipublic attributes on support contingent on familiarity and experiences - US sample
diff_amce_exp_choice_usa <- cj(experiences_usa, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~familiar)
diff_amce_participated_choice_usa <- cj(experiences_usa, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~participated)
diff_amce_expdirect_choice_usa <- cj(usa_expdirect, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~exsplit)

diff_experience_estimates_choice_usa <-bind_rows(
  diff_amce_exp_choice_usa %>% mutate(model = "Familiarity"),
  diff_amce_participated_choice_usa %>% mutate(model = "Participated"),
  diff_amce_expdirect_choice_usa %>% mutate(model = "Positive experience"))

diff_experience_choice_usa <- diff_experience_estimates_choice_usa %>% select(BY:level, estimate:model)
diff_experience_choice_usa$model <- factor(diff_experience_choice_usa$model, levels = c("Familiarity", "Participated", "Positive experience"))

plot(diff_experience_choice_usa, feature_headers = FALSE, xlab = "Estimated Difference (AMCE)") + 
  ggplot2::facet_wrap(~model, ncol = 3) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("USA") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 17: AMCE of minipublic attributes on support contingent on familiarity and experiences - US sample
head(diff_experience_choice_usa[c("level", "estimate", "std.error")], 30L)

# Figure D. 17: MARGINAL MEANS for familiarity and experiences – US sample
diff_mm_exp_choice_usa <- cj(experiences_usa, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~familiar)
diff_mm_participated_choice_usa <- cj(experiences_usa, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~participated)
diff_mm_expdirect_choice_usa <- cj(usa_expdirect, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~exsplit)

diff_experience_estimates_choice_usa_mm <-bind_rows(
  diff_mm_exp_choice_usa %>% mutate(model = "Familiarity"),
  diff_mm_participated_choice_usa %>% mutate(model = "Participated"),
  diff_mm_expdirect_choice_usa %>% mutate(model = "Positive experience"))

diff_experience_choice_usa_mm <- diff_experience_estimates_choice_usa_mm %>% select(BY:level, estimate:model)
diff_experience_choice_usa_mm$model <- factor(diff_experience_choice_usa_mm$model, levels = c("Familiarity", "Participated", "Positive experience"))

plot(diff_experience_choice_usa_mm, feature_headers = FALSE, xlab = "Estimated Difference (Marginal Means)") + 
  ggplot2::facet_wrap(~model, ncol = 3) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("USA") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 18: MARGINAL MEANS for familiarity and experiences – US sample
head(diff_experience_choice_usa_mm[c("level", "estimate", "std.error")], 57L)

#### Trust ####
## Reproduce: 
##      - Figure D.18 (Appendix)
##      - Table D.19 (Appendix)
##      - Figure D.19 (Appendix)
##      - Table D.20 (Appendix)
##
## Output: 
##      Figure D.18: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citi-zens (bottom row) – US sample
##	    Table D.19: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citizens (bottom row) – US sample
##	    Figure D.19: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – US sample
##	    Table D.20: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – US sample
##
# Trust in political institutions: index parliament and government
alpha(subset(cj_usa, select = c(trust_parliament, trust_government)), check.keys = TRUE) # alpha .91
cj_usa$trust_politics <- rowMeans(subset(cj_usa, select = c(trust_parliament, trust_government)), na.rm = TRUE)
describe(cj_usa$trust_politics)
cj_usa_politics <- cj_usa %>% filter(trust_politics != "NaN")
describe(cj_usa$trust_politics) #mean 3.36, median 3
cj_usa_politics$trustsplit<-dicho(cj_usa_politics$trust_politics, dich.by = "median", as.num = TRUE)
class(cj_usa_politics$trustsplit)
cj_usa_politics$trustsplit <- factor(cj_usa_politics$trustsplit, labels = c("Low trust", "High trust"))
cj_usa_politics$trustsplit <- relevel(cj_usa_politics$trustsplit, ref = "High trust")

summary(cj_usa_politics$trustsplit)

cj_anova(cj_usa_politics, Mchoice, by = ~trustsplit)

# Trust in politicians
usa_politicians_NA <- is.na(select(cj_usa, trust_politicians))
cj_usa_politicians <- filter(cj_usa, usa_politicians_NA == 0)
describe(cj_usa_politicians$trust_politicians) #mean 2.86 median 3
cj_usa_politicians$trustsplit<-dicho(cj_usa_politicians$trust_politicians, dich.by = "median", as.num = TRUE)
class(cj_usa_politicians$trustsplit)
cj_usa_politicians$trustsplit <- factor(cj_usa_politicians$trustsplit, labels = c("Low trust", "High trust"))
cj_usa_politicians$trustsplit <- relevel(cj_usa_politicians$trustsplit, ref = "High trust")

summary(cj_usa_politicians$trustsplit)

cj_anova(cj_usa_politicians, Mchoice, by = ~trustsplit)

# Trust in minipublics 
describe(cj_usa$trust_dmp)#mean: 4.58 median: 5
cj_usa$dmpsplit<-dicho(cj_usa$trust_dmp, dich.by = "median", as.num = TRUE)
class(cj_usa$dmpsplit)
cj_usa$dmpsplit <- factor(cj_usa$dmpsplit, labels = c("Low trust", "High trust"))

summary(cj_usa$dmpsplit)

cj_anova(cj_usa, Mchoice, by = ~dmpsplit)

# Trust in other citizens
alpha(subset(cj_usa, select = c(trust_citizens1, trust_citizens2)), check.keys = TRUE) # alpha .87
cj_usa$trust_citizens <- rowMeans(subset(cj_usa, select = c(trust_citizens1, trust_citizens2)), na.rm = TRUE)
describe(cj_usa$trust_citizens) #mean 4.42 median 4.5
cj_usa$citizensplit<-dicho(cj_usa$trust_citizens, dich.by = "median", as.num = TRUE)
class(cj_usa$citizensplit)
cj_usa$citizensplit <- factor(cj_usa$citizensplit, labels = c("Low trust", "High trust"))

summary(cj_usa$citizensplit)

cj_anova(cj_usa, Mchoice, by = ~citizensplit)

describeBy(cj_usa_politics$retro, cj_usa_politics$trustsplit)
describeBy(cj_usa_politicians$retro, cj_usa_politicians$trustsplit)
describeBy(cj_usa$retro, cj_usa$dmpsplit)
describeBy(cj_usa$retro, cj_usa$citizensplit)

t.test(cj_usa_politics$retro ~ cj_usa_politics$trustsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.28; low trust: mean=4.81)
cohensD(cj_usa_politics$retro ~ cj_usa_politics$trustsplit) # Cohens d = 0.34 (small effect)

t.test(cj_usa_politicians$retro ~ cj_usa_politicians$trustsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.40; low trust: mean=4.86)
cohensD(cj_usa_politicians$retro ~ cj_usa_politicians$trustsplit) # Cohens d = 0.39 (small effect)

t.test(cj_usa$retro ~ cj_usa$dmpsplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=6.28; low trust: mean=4.63)
cohensD(cj_usa$retro ~ cj_usa$dmpsplit) # Cohens d = 1.39 (large effect)

t.test(cj_usa$retro ~ cj_usa$citizensplit, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.80; low trust: mean=4.46)
cohensD(cj_usa$retro ~ cj_usa$citizensplit) # Cohens d = 1.11 (large effect)

# Figure D. 18: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citi-zens (bottom row) – US sample
usa_diff_amce_politics <- cj(cj_usa_politics, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
usa_diff_amce_politicians <- cj(cj_usa_politicians, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
usa_diff_amce_dmp <- cj(cj_usa, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~dmpsplit)
usa_diff_amce_citizens <- cj(cj_usa, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~citizensplit)

usa_diff_disent_estimates_trust <-bind_rows(
  usa_diff_amce_politics %>% mutate(model = "Political institutions"),
  usa_diff_amce_politicians %>% mutate(model = "Politicians"),
  usa_diff_amce_dmp %>% mutate(model = "Minipublics"),
  usa_diff_amce_citizens %>% mutate(model = "Citizens"))

usa_diff_disent_trust <- usa_diff_disent_estimates_trust %>% select(BY:level, estimate:model)
usa_diff_disent_trust$model <- factor(usa_diff_disent_trust$model, levels = c("Political institutions", "Politicians", "Minipublics", "Citizens"))

plot(usa_diff_disent_trust, feature_headers = FALSE, xlab = "Estimated Difference (AMCE)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("USA") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 19: AMCE of minipublic attributes on support contingent on low political trust (top row) and high trust in citizens (bottom row) – US sample
head(usa_diff_disent_trust[c("level", "estimate", "std.error")], 40L)

# Figure D. 19: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – US sample
usa_diff_mm_politics <- cj(cj_usa_politics, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
usa_diff_mm_politicians <- cj(cj_usa_politicians, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit)
usa_diff_mm_dmp <- cj(cj_usa, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~dmpsplit)
usa_diff_mm_citizens <- cj(cj_usa, Mchoice, id = ~ID, estimate="mm_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~citizensplit)

usa_diff_disent_estimates_trust_mm <-bind_rows(
  usa_diff_mm_politics %>% mutate(model = "Political institutions"),
  usa_diff_mm_politicians %>% mutate(model = "Politicians"),
  usa_diff_mm_dmp %>% mutate(model = "Minipublics"),
  usa_diff_mm_citizens %>% mutate(model = "Citizens"))

usa_diff_disent_trust_mm <- usa_diff_disent_estimates_trust_mm %>% select(BY:level, estimate:model)
usa_diff_disent_trust_mm$model <- factor(usa_diff_disent_trust_mm$model, levels = c("Political institutions", "Politicians", "Minipublics", "Citizens"))

plot(usa_diff_disent_trust_mm, feature_headers = FALSE, xlab = "Estimated Difference (Marginal Means)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  ggtitle("USA") +
  scale_color_manual(values = rep("#515151", 9))

# Table D. 20: MARGINAL MEANS for low political trust (top row) and high trust in citizens (bottom row) – US sample
head(usa_diff_disent_trust_mm[c("level", "estimate", "std.error")], 76L)

######### APPENDIX E: Further robustness and attention checks #########
#### Rating outcome variable ####
# (note: no data available for Finland)
## Reproduce: 
##      - Figure E.1 (Appendix)
##      - Table E.1 (Appendix)
##
## Output: 
##      Figure E.1: Effects of minipublic attributes on support - rating outcome variable
##	    Table E.1: Effects of minipublic attributes on support - rating outcome variable
##
describe(cj_all$rating)
Mrating <- rating ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav # conjoint with all attributes, outcome favorability, rating outcome
amce_rating <- cregg::cj(cj_all, Mrating, id = ~ID,  estimate="amce", feature_labels = feature_names)

# Figure E. 1: Effects of minipublic attributes on support - rating outcome variable
plot(amce_rating, feature_headers = FALSE)+
  ggplot2::facet_grid(feature~.,  scales = "free_y", space = "free", switch ="y") +
  theme(panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 12),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none",
        panel.spacing.y = unit(1, "lines"),
        strip.background= element_rect(color="white", size = 2, fill = "gray97"), 
        strip.text.y.left = element_text(angle = 0),
        strip.placement = "outside",
        strip.text.y = element_text(colour = "#515151", size = 10, face = "bold")) +
  theme(plot.title = element_text(colour =  "#515151", size = 14)) +
  geom_point(size=2) +
  scale_color_manual(values = rep("#515151", 9))

# Table E. 1: Effects of minipublic attributes on support - rating outcome variable
C20 <- head(amce_rating [c("level", "estimate", "std.error")], 19L)
knitr::kable(C20, format = "markdown")

#### Marginal Means plot ####
## Reproduce: 
##      - Figure E.2 (Appendix)
##      - Table E.2 (Appendix)
##
## Output: 
##      Figure E.2: Effects of minipublic attributes on support – Marginal means
##	    Table E.2: Effects of minipublic attributes on support – Marginal means
##
mm_choice <- cj(cj_all, Mchoice, id = ~ID, estimate="mm", h0=0.5, feature_labels = feature_names)

# Figure E. 2: Effects of minipublic attributes on support – Marginal means 
plot(mm_choice, feature_headers = FALSE)+
  ggplot2::facet_grid(feature~.,  scales = "free_y", space = "free", switch ="y") +
  theme(panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 12),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none",
        panel.spacing.y = unit(1, "lines"),
        strip.background= element_rect(color="white", size = 2, fill = "gray97"), 
        strip.text.y.left = element_text(angle = 0),
        strip.placement = "outside",
        strip.text.y = element_text(colour = "#515151", size = 10, face = "bold")) +
  theme(plot.title = element_text(colour =  "#515151", size = 14)) +
  geom_point(size=2) +
  geom_vline(xintercept = 0.5, colour="grey") +
  scale_color_manual(values = rep("#515151", 9))

# Table E.2: Effects of minipublic attributes on support – Marginal means 
C21 <- head(mm_choice [c("level", "estimate", "std.error")], 19L)
knitr::kable(C21, format = "markdown")

#### Satisfaction w/ democracy ####
## Reproduce: 
##      - Figure E.3 (Appendix)
##      - Table E.3 (Appendix)
##
## Output: 
##      Figure E.3: AMCE for dissatisfied and satisfied citizens
##	    Table E.3: AMCE for dissatisfied and satisfied citizens
##
describe(cj_trust$swd)
cj_trust$swd[cj_trust$swd==997] <- NA
cj_trust$swd[cj_trust$swd==8] <- NA
describe(cj_trust$swd) # mean 4.41 median 5

# prepare grouping variable for swd:
# median split (use dichotomize function)
cj_trust$swdsplit<-dicho(cj_trust$swd, dich.by = "median", as.num = TRUE)
cj_trust$swdsplit <- factor(cj_trust$swdsplit, labels = c("dissatisfaction w/ democracy", "satisfaction w/ democracy"))
cj_trust$swdsplit <- relevel(cj_trust$swdsplit, ref = "satisfaction w/ democracy")

# (!) use different cutoffs: use median split rule used in Goldberg and Bächtiger, 2023
# If median > mean then ‘low’ is assigned for < median and ‘high’ for ≥ median; 
# If median < mean then ‘low’ is assigned for ≤ median and ‘high’ for > median.

cj_trust$swdsplit_2[cj_trust$swd >= 5] <- 1L
cj_trust$swdsplit_2[cj_trust$swd < 5] <- 2L
cj_trust$swdsplit_2 <- factor(cj_trust$swdsplit_2, 1:2, c("satisfied", "dissatisfied"))
summary(cj_trust$swdsplit_2)

# Figure E. 3: AMCE for dissatisfied and satisfied citizens
amce_by_sat_choice <- cj(cj_trust, Mchoice, id = ~ID,estimate = "amce",feature_labels = feature_names, by = ~swdsplit_2)
diff_amce_sat_cond_choice <- cj(cj_trust, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, by = ~swdsplit_2)
amce_sat_choice_diff <- rbind(amce_by_sat_choice, diff_amce_sat_cond_choice)

plot(amce_sat_choice_diff, feature_headers = FALSE) +
  ggplot2::facet_grid(feature~BY, scales = "free_y", space="free", switch = "y") +
  theme(strip.background.x = element_blank(),
        strip.text.x = element_text(colour = "#091926", size = 14),
        strip.text.y = element_blank(),
        panel.spacing.y = unit(1, "lines"),
        panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#091926", size = 12),
        axis.title = element_text(colour = "#091926", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = c("#515151", "#515151"))

# Table E. 3: AMCE for dissatisfied and satisfied citizens
C22 <- head(amce_sat_choice_diff [c("level", "estimate", "std.error")], 48L)
knitr::kable(C22, format = "markdown")

#### Different cutoffs ####
# use different cutoffs: use median split rule used in Goldberg and Bächtiger, 2023
# If median > mean then ‘low’ is assigned for < median and ‘high’ for ≥ median; 
# If median < mean then ‘low’ is assigned for ≤ median and ‘high’ for > median.
## Reproduce: 
##      - Figure E.4 (Appendix)
##      - Table E.4 (Appendix)
##      - Figure E.5 (Appendix)
##      - Table E.5 (Appendix)
##      - Figure E.6 (Appendix)
##      - Table E.6 (Appendix)
##
## Output: 
##      Figure E.4: Effects of minipublic attributes on support conditional on low political trust (top row) and high trust in other citizens (bot-tom row) – different cutoffs
##	    Table E.4: Effects of minipublic attributes on support conditional on low political trust (top row) and high trust in other citizens (bottom row) – different cutoffs
##	    Figure E.5: Combined effect for low political trust but high trust in citizens
##      Table E.5: Combined effect for low political trust but high trust in citizens
##      Figure E.6: Effects of minipublic attributes on support conditional on familiarity and positive experiences – different cutoffs
##      Table E.6: Effects of minipublic attributes on support conditional on familiarity and positive experiences – different cutoffs
##
# Trust in politics
cj_trust_politics$trustsplit_b[cj_trust_politics$trust_politics >= 4] <- 1L
cj_trust_politics$trustsplit_b[cj_trust_politics$trust_politics < 4] <- 2L
cj_trust_politics$trustsplit_b <- factor(cj_trust_politics$trustsplit_b, 1:2, c("High trust", "Low trust"))
# nested model comparison test (does any of the interactions between the "by variable" and feature levels differ from zero?)
cj_anova(cj_trust_politics, Mchoice, by = ~trustsplit_b)

summary(cj_trust_politics$trustsplit_b)

# Trust in politicians 
cj_trust_politicians$trustsplit_b[cj_trust_politicians$trust_politicians > 3] <- 1L
cj_trust_politicians$trustsplit_b[cj_trust_politicians$trust_politicians <= 3] <- 2L
cj_trust_politicians$trustsplit_b <- factor(cj_trust_politicians$trustsplit_b, 1:2, c("High trust", "Low trust"))
# nested model comparison test (does any of the interactions between the "by variable" and feature levels differ from zero?)
cj_anova(cj_trust_politicians, Mchoice, by = ~trustsplit_b)

summary(cj_trust_politicians$trustsplit_b)

# Trust in minipublics
cj_trust$dmpsplit_b[cj_trust$trust_dmp >= 5] <- 1L
cj_trust$dmpsplit_b[cj_trust$trust_dmp < 5] <- 2L
cj_trust$dmpsplit_b <- factor(cj_trust$dmpsplit_b, 1:2, c("High trust", "Low trust"))
cj_trust$dmpsplit_b <- relevel(cj_trust$dmpsplit_b, ref = "Low trust")
# nested model comparison test (does any of the interactions between the "by variable" and feature levels differ from zero?)
cj_anova(cj_trust, Mchoice, by = ~dmpsplit_b)

summary(cj_trust$dmpsplit_b)

# Trust in other citizens
cj_trust$citizensplit_b[cj_trust$trust_citizens >= 4.5] <- 1L
cj_trust$citizensplit_b[cj_trust$trust_citizens < 4.5] <- 2L
cj_trust$citizensplit_b <- factor(cj_trust$citizensplit_b, 1:2, c("High trust", "Low trust"))
cj_trust$citizensplit_b <- relevel(cj_trust$citizensplit_b, ref = "Low trust")
# nested model comparison test (does any of the interactions between the "by variable" and feature levels differ from zero?)
cj_anova(cj_trust, Mchoice, by = ~citizensplit_b)

summary(cj_trust$citizensplit_b)

# Figure E. 4: Effects of minipublic attributes on support conditional on low political trust (top row) and high trust in other citizens (bottom row) – different cutoffs
diff_amce_politics_b <- cj(cj_trust_politics, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit_b)
diff_amce_politicians_b <- cj(cj_trust_politicians, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~trustsplit_b)
diff_amce_dmp_b <- cj(cj_trust, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~dmpsplit_b)
diff_amce_citizens_b <- cj(cj_trust, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~citizensplit_b)

diff_disent_estimates_trust_b <-bind_rows(
  diff_amce_politics_b %>% mutate(model = "Political institutions"),
  diff_amce_politicians_b %>% mutate(model = "Politicians"),
  diff_amce_dmp_b %>% mutate(model = "Minipublics"),
  diff_amce_citizens_b %>% mutate(model = "Citizens"))

diff_disent_trust_b <- diff_disent_estimates_trust_b %>% select(BY:level, estimate:model)
diff_disent_trust_b$model <- factor(diff_disent_trust_b$model, levels = c("Political institutions", "Politicians", "Minipublics", "Citizens"))

plot(diff_disent_trust_b, feature_headers = FALSE, xlab = "Estimated Difference (AMCE)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = rep("#515151", 9))

# Table E. 4: Effects of minipublic attributes on support conditional on low political trust (top row) and high trust in other citizens (bottom row) – different cutoffs
C23 <- head(diff_disent_trust_b [c("level", "estimate", "std.error")], 40L)
knitr::kable(C23, format = "markdown")

t.test(cj_trust_politics$retro ~ cj_trust_politics$trustsplit_b, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.00; low trust: mean=4.69)
cohensD(cj_trust_politics$retro ~ cj_trust_politics$trustsplit_b) # Cohens d = 0.23 (small effect)

t.test(cj_trust_politicians$retro ~ cj_trust_politicians$trustsplit_b, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5; low trust: mean=4.74)
cohensD(cj_trust_politicians$retro ~ cj_trust_politicians$trustsplit_b) # Cohens d = 0.19 (small effect)

t.test(cj_trust$retro ~ cj_trust$dmpsplit_b, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.58; low trust: mean=4.05)
cohensD(cj_trust$retro ~ cj_trust$dmpsplit_b) # Cohens d = 1.35 (large effect)

t.test(cj_trust$retro ~ cj_trust$citizensplit_b, var.equal = TRUE, alternative = "two.sided") # significant differences (high trust: mean=5.42; low trust: mean=4.17)
cohensD(cj_trust$retro ~ cj_trust$citizensplit_b) # Cohens d = 1.04 (large effect)

# Figure E. 5: Combined effect for low political trust but high trust in citizens
# all countries
# split scales at midpoint (<= and >=) to define cases where political trust is low but trust in citizens is high 
cj_trust_NA$select_trust_b <- ifelse(cj_trust_NA$poltrust <= 4 & cj_trust_NA$soctrust >= 4, 1, 0)
prop.table(table(cj_trust_NA$select_trust_b)) #57,6% in group 0 and 42,4% in group 1
cj_trust_NA$select_trust_b <- as.factor(cj_trust_NA$select_trust_b)
summary(cj_trust_NA$select_trust_b)
cj_anova(cj_trust_NA, Mchoice, by = ~select_trust_b)

# Finland
# split scales at midpoint to define cases where political trust is low but trust in citizens is high 
cj_fin_NA$select_trust_b <- ifelse(cj_fin_NA$poltrust <= 4 & cj_fin_NA$soctrust >= 4, 1, 0)
prop.table(table(cj_fin_NA$select_trust_b)) #67,2% in group 0 and 32,8% in group 1
cj_fin_NA$select_trust_b <- as.factor(cj_fin_NA$select_trust_b)
summary(cj_fin_NA$select_trust_b)
cj_anova(cj_fin_NA, Mchoice, by = ~select_trust_b)

# Ireland
# split scales at midpoint to define cases where political trust is low but trust in citizens is high 
cj_ire_NA$select_trust_b <- ifelse(cj_ire_NA$poltrust <= 4 & cj_ire_NA$soctrust >= 4, 1, 0)
prop.table(table(cj_ire_NA$select_trust_b)) #55,9% in group 0 and 44,1% in group 1
cj_ire_NA$select_trust_b <- as.factor(cj_ire_NA$select_trust_b)
summary(cj_ire_NA$select_trust_b)
cj_anova(cj_ire_NA, Mchoice, by = ~select_trust_b)

# USA
# split scales at midpoint to define cases where political trust is low but trust in citizens is high 
cj_usa_NA$select_trust_b <- ifelse(cj_usa_NA$poltrust <= 4 & cj_usa_NA$soctrust >= 4, 1, 0)
prop.table(table(cj_usa_NA$select_trust_b)) #50% in group 0 and 50% in group 1
cj_usa_NA$select_trust_b <- as.factor(cj_usa_NA$select_trust_b)
summary(cj_usa_NA$select_trust_b)
cj_anova(cj_usa_NA, Mchoice, by = ~select_trust_b)

diff_amce_select_all_b <- cj(cj_trust_NA, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~select_trust_b)
diff_amce_select_fin_b <- cj(cj_fin_NA, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~select_trust_b)
diff_amce_select_ire_b <- cj(cj_ire_NA, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~select_trust_b)
diff_amce_select_usa_b <- cj(cj_usa_NA, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~select_trust_b)

diff_disent_estimates_select_b <-bind_rows(
  diff_amce_select_all_b %>% mutate(model = "All countries"),
  diff_amce_select_fin_b %>% mutate(model = "Finland"),
  diff_amce_select_ire_b %>% mutate(model = "Ireland"),
  diff_amce_select_usa_b %>% mutate(model = "USA"))

diff_disent_select_b <- diff_disent_estimates_select_b %>% select(BY:level, estimate:model)
diff_disent_select_b$model <- factor(diff_disent_select_b$model, levels = c("All countries", "Finland", "Ireland", "USA"))

plot(diff_disent_select_b, feature_headers = FALSE, xlab = "Estimated Difference (AMCE)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background= element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = rep("#515151", 9))

#Table E. 5: Combined effect for low political trust but high trust in citizens
C24 <- head(diff_disent_select_b [c("level", "estimate", "std.error")], 40L)
knitr::kable(C24, format = "markdown")

# Experience
describe(cj_expdirect$experiences_direct) # mean 4.72 median 5
cj_expdirect$exsplit_b[cj_expdirect$experiences_direct >= 5] <- 1L
cj_expdirect$exsplit_b[cj_expdirect$experiences_direct < 5] <- 2L
cj_expdirect$exsplit_b <- factor(cj_expdirect$exsplit_b, 1:2, c("Positive experiences", "Negative experiences"))
cj_expdirect$exsplit_b <- relevel(cj_expdirect$exsplit_b, ref = "Negative experiences")
# nested model comparison test (does any of the interactions between the "by variable" and feature levels differ from zero?)
cj_anova(cj_expdirect, Mchoice, by = ~exsplit_b)

summary(cj_expdirect$exsplit_b)

# Figure E. 6: Effects of minipublic attributes on support conditional on familiarity and positive experiences – different cutoffs
diff_amce_exp_choice_b <- cj(cj_experience, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~familiar)
diff_amce_expdirect_choice_b <- cj(cj_expdirect, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~exsplit_b)

diff_experience_estimates_choice_b <-bind_rows(
  diff_amce_exp_choice_b %>% mutate(model = "Familiarity"),
  diff_amce_expdirect_choice_b %>% mutate(model = "Positive experience"))

diff_experience_choice_b <- diff_experience_estimates_choice_b %>% select(BY:level, estimate:model)
diff_experience_choice_b$model <- factor(diff_experience_choice_b$model, levels = c("Familiarity", "Positive experience"))

plot(diff_experience_choice_b, feature_headers = FALSE, xlab = "Estimated Difference (AMCE)") + 
  ggplot2::facet_wrap(~model, ncol = 2) +
  theme(panel.border = element_rect(linetype = "dotted", colour = "#515151"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 10),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        strip.background = element_rect(color="gray97", size = 2, fill = "gray97"),
        strip.text.x = element_text(colour = "#515151", size = 14),
        legend.position = "none")+
  geom_point(size=2) +
  scale_color_manual(values = rep("#515151", 9))

# Table E. 6: Effects of minipublic attributes on support conditional on familiarity and positive experiences – different cutoffs
C25 <- head(diff_experience_choice_b [c("level", "estimate", "std.error")], 20L)
knitr::kable(C25, format = "markdown")

#### Attention checks ####
## Reproduce: 
##      - Figure E.7 (Appendix)
##      - Table E.7 (Appendix)
##      - Figure E.8 (Appendix)
##      - Table E.8 (Appendix)
##      - Figure E.9 (Appendix)
##      - Table E.9 (Appendix)
##
## Output: 
##      Figure E.7: Attention check – video watched
##	    Table E.7: Attention check – video watched
##	    Figure E.8: Attention check – reading arguments (authorization)
##      Table E.8: Attention check – reading arguments (authorization)
##      Figure E.9: Attention check – reading arguments (recruitment)
##      Table E.9: Attention check – reading arguments (recruitment)
##
# relevant variables: time spent watching the video (all countries), time spent reading the arguments (IRE and US), time spent for the conjoint tasks (IRE and US)
describe(cj_trust$video_time_spent)
# drop outliers (video was about 90 seconds)
cj_video_watched <- cj_trust |>  dplyr::filter(video_time_spent <= 120)
describe(cj_video_watched$video_time_spent) # mean 60.17 median 90
# prepare grouping variable at halftime of video length
cj_video_watched$split[cj_video_watched$video_time_spent > 45] <- 1L
cj_video_watched$split[cj_video_watched$video_time_spent <= 45] <- 2L
cj_video_watched$split <- factor(cj_video_watched$split, 1:2, c("at least half of the video was watched", "less than half od the video was watched"))
# nested model comparison test (does any of the interactions between the "by variable" and feature levels differ from zero?)
cj_anova(cj_video_watched, Mchoice, by = ~split)
summary(cj_video_watched$split)

# Figure E. 7: Attention check – video watched
amce_by_video <- cj(cj_video_watched, Mchoice, id = ~ID, estimate = "amce", feature_labels = feature_names, feature_order = c("outcome.fav",  "authorization", "aim", "majority", "composition",  "size", "recruitment", "initiative", "issue"), by = ~split)

plot(amce_by_video, group = "split", feature_headers = FALSE)+
  ggplot2::facet_wrap(~BY, ncol = 2L) +
  theme(panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 12),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none",
        panel.spacing.y = unit(1, "lines"),
        strip.background= element_rect(color="white", size = 2, fill = "gray97"), 
        strip.text.y.left = element_text(angle = 0),
        strip.placement = "outside",
        strip.text.y = element_text(colour = "#515151", size = 10, face = "bold")) +
  theme(plot.title = element_text(colour =  "#515151", size = 14)) +
  geom_point(size=2) +
  scale_color_manual(values = c("#515151", "#515151", "#515151"))

# Table E. 7: Attention check – video watched
C26 <- head(amce_by_video [c("level", "estimate", "std.error")], 38L)
knitr::kable(C26, format = "markdown")

# reading arguments (note: data only available for Ireland and Finland)
cj_attention_IREUS <- cj_trust |>  dplyr::filter(country != "FIN")

describe(cj_attention_IREUS$recruitment_time_spent)
describe(cj_attention_IREUS$authorization_time_spent)

# identify Q1, Q3 and IQR
#The interquartile range (IQR) is the difference between the 75th percentile (Q3) 
#and the 25th percentile (Q1) in a data set. It measures the dispersion of the middle 50% of the values.
Q1_A <- quantile(cj_attention_IREUS$authorization_time_spent, .25)
Q3_A <- quantile(cj_attention_IREUS$authorization_time_spent, .75)
IQR_A <- IQR(cj_attention_IREUS$authorization_time_spent)

Q1_R <- quantile(cj_attention_IREUS$recruitment_time_spent, .25)
Q3_R <- quantile(cj_attention_IREUS$recruitment_time_spent, .75)
IQR_R <- IQR(cj_attention_IREUS$recruitment_time_spent)

# keep rows with values which are 1,5 * IQR in Q1 and Q3 
no_outliers_A <- subset(cj_attention_IREUS, cj_attention_IREUS$authorization_time_spent > (Q1_A - 1.5*IQR_A) & cj_attention_IREUS$authorization_time_spent < (Q3_A + 1.5*IQR_A))
no_outliers_R <- subset(cj_attention_IREUS, cj_attention_IREUS$recruitment_time_spent > (Q1_R - 1.5*IQR_R) & cj_attention_IREUS$recruitment_time_spent < (Q3_R + 1.5*IQR_R))

summary(no_outliers_A$authorization_time_spent)
describe(no_outliers_A$authorization_time_spent)
summary(no_outliers_R$recruitment_time_spent)
describe(no_outliers_R$recruitment_time_spent)

# read authorization (split at mean)
no_outliers_A$readAsplit[no_outliers_A$authorization_time_spent < 12.81] <- 1L
no_outliers_A$readAsplit[no_outliers_A$authorization_time_spent >= 12.81] <- 2L
no_outliers_A$readAsplit <- factor(no_outliers_A$readAsplit, 1:2, c("scarcely read", "thoroughly read"))
summary(no_outliers_A$readAsplit)
cj_anova(no_outliers_A, Mchoice, by = ~readAsplit) # sig

# read recruitment (split at mean)
no_outliers_R$readRsplit[no_outliers_R$recruitment_time_spent < 14.36] <- 1L
no_outliers_R$readRsplit[no_outliers_R$recruitment_time_spent >= 12.36] <- 2L
no_outliers_R$readRsplit <- factor(no_outliers_R$readRsplit, 1:2, c("scarcely read", "thoroughly read"))
summary(no_outliers_R$readRsplit)
cj_anova(no_outliers_R, Mchoice, by = ~readRsplit) # sig

# read authorization (split at median)
no_outliers_A$readAsplit_b[no_outliers_A$authorization_time_spent < 6] <- 1L
no_outliers_A$readAsplit_b[no_outliers_A$authorization_time_spent >= 6] <- 2L
no_outliers_A$readAsplit_b <- factor(no_outliers_A$readAsplit_b, 1:2, c("scarcely read", "thoroughly read"))
summary(no_outliers_A$readAsplit_b)
cj_anova(no_outliers_A, Mchoice, by = ~readAsplit_b) # sig

# read recruitment (split at median)
no_outliers_R$readRsplit_b[no_outliers_R$recruitment_time_spent < 9] <- 1L
no_outliers_R$readRsplit_b[no_outliers_R$recruitment_time_spent >= 9] <- 2L
no_outliers_R$readRsplit_b <- factor(no_outliers_R$readRsplit_b, 1:2, c("scarcely read", "thoroughly read"))
summary(no_outliers_R$readRsplit_b)
cj_anova(no_outliers_R, Mchoice, by = ~readRsplit_b) # sig

# Figure E. 8: Attention check – reading arguments (authorization)
amce_by_authorization <- cj(no_outliers_A, Mchoice, id = ~ID, estimate = "amce", feature_labels = feature_names, feature_order = c("outcome.fav",  "authorization", "aim", "majority", "composition",  "size", "recruitment", "initiative", "issue"), by = ~readAsplit)

plot(amce_by_authorization, group = "readAsplit", feature_headers = FALSE)+
  ggplot2::facet_wrap(~BY, ncol = 2L) +
  theme(panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 12),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none",
        panel.spacing.y = unit(1, "lines"),
        strip.background= element_rect(color="white", size = 2, fill = "gray97"), 
        strip.text.y.left = element_text(angle = 0),
        strip.placement = "outside",
        strip.text.y = element_text(colour = "#515151", size = 10, face = "bold")) +
  theme(plot.title = element_text(colour =  "#515151", size = 14)) +
  geom_point(size=2) +
  scale_color_manual(values = c("#515151", "#515151", "#515151"))

# Table E. 8: Attention check – reading arguments (authorization)
C27 <- head(amce_by_authorization [c("level", "estimate", "std.error")], 38L)
knitr::kable(C27, format = "markdown")

# Figure E. 9: Attention check – reading arguments (recruitment)
amce_by_recruitment <- cj(no_outliers_R, Mchoice, id = ~ID, estimate = "amce", feature_labels = feature_names, feature_order = c("outcome.fav",  "authorization", "aim", "majority", "composition",  "size", "recruitment", "initiative", "issue"), by = ~readRsplit)

plot(amce_by_recruitment, group = "readRsplit", feature_headers = FALSE)+
  ggplot2::facet_wrap(~BY, ncol = 2L) +
  theme(panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#515151", size = 12),
        axis.title = element_text(colour = "#515151", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none",
        panel.spacing.y = unit(1, "lines"),
        strip.background= element_rect(color="white", size = 2, fill = "gray97"), 
        strip.text.y.left = element_text(angle = 0),
        strip.placement = "outside",
        strip.text.y = element_text(colour = "#515151", size = 10, face = "bold")) +
  theme(plot.title = element_text(colour =  "#515151", size = 14)) +
  geom_point(size=2) +
  scale_color_manual(values = c("#515151", "#515151", "#515151"))

# Table E. 9: Attention check – reading arguments (recruitment)
C28 <- head(amce_by_recruitment [c("level", "estimate", "std.error")], 38L)
knitr::kable(C28, format = "markdown")

#### Carryover, balancing, and profile order effects ####
## Reproduce: 
##      - Figure E.10 (Appendix)
##      - Figure E.11 (Appendix)
##      - Table E.10 (Appendix)
##      - Table E.11 (Appendix)
##      - Figure E.12 (Appendix)
##      - Table E.12 (Appendix)
##
## Output: 
##      Figure E.10: Test for carry-over effects
##	    Figure E.11: Test for randomization / balancing
##	    Table E.10: Test for carry-over effects
##      Table E.11: Test for randomization / balancing
##      Figure E.12: Test for profile order effects
##      Table E.12: Test for profile order effects
##
cj_all$task <- factor(cj_all$task)

amce_by_tasks <- cj(cj_all, Mchoice, id = ~ID, estimate = "amce", feature_labels = feature_names, feature_order = c("outcome.fav",  "authorization", "aim", "majority", "composition",  "size", "recruitment", "initiative", "issue"), by = ~task)
summary(cj_all$task)

tasks <- plot(cj(cj_all, Mchoice, id = ~ID, by = ~task, estimate = "amce"), group = "task")

# Figure E. 10: Test for carry-over effects
c <- tasks + theme(panel.border = element_rect(linetype = "blank"),
                   panel.background = element_blank(), 
                   axis.text.y = element_text(colour = "#515151", size = 12),
                   axis.text.x = element_text(colour = "#091926", size = 12),
                   axis.title = element_text(colour = "#091926", size = 12),
                   panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
                   axis.ticks = element_blank(), 
                   legend.position = "none",
                   panel.spacing.y = unit(1, "lines"),
                   strip.background= element_rect(color="white", size = 2, fill = "#f9f9f9"), 
                   strip.text.y.left = element_text(angle = 0),
                   strip.placement = "outside",
                   strip.text.y = element_text(colour = "#091926", size = 12, face = "bold")) +
  xlab("Change: Pr(choice DMP)") +
  scale_color_manual(values = rep("#515151", 9))

# only first task 
cj_tasks_1 <- cj_all %>%  filter(task == "1")
choice_amce_tasks_1 <- amce(cj_tasks_1, Mchoice, id = ~ID, feature_labels = feature_names)

c1 <- plot(choice_amce_tasks_1, feature_headers = FALSE)+
  ggplot2::facet_grid(feature~.,  scales = "free_y", space = "free", switch ="y") +
  theme(panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#091926", size = 12),
        axis.title = element_text(colour = "#091926", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none",
        panel.spacing.y = unit(1, "lines"),
        strip.background= element_rect(color="white", size = 2, fill = "#f9f9f9"), 
        strip.text.y.left = element_text(angle = 0),
        strip.placement = "outside",
        strip.text.y = element_text(colour = "#091926", size = 12, face = "bold")) +
  xlab("Change: Pr(choice DMP)") +
  ggtitle("Task 1") +
  geom_point(size=2) +
  scale_color_manual(values = rep("#515151", 9))

# exclude first task
cj_tasks <- cj_all %>%  filter(task != "1")
choice_amce_tasks <- amce(cj_tasks, Mchoice, id = ~ID, feature_labels = feature_names)
rating_amce_tasks <- amce(cj_tasks, Mrating, id = ~ID, feature_labels = feature_names)

c2 <- plot(choice_amce_tasks, feature_headers = FALSE)+
  ggplot2::facet_grid(feature~.,  scales = "free_y", space = "free", switch ="y") +
  theme(panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#091926", size = 12),
        axis.title = element_text(colour = "#091926", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none",
        panel.spacing.y = unit(1, "lines"),
        strip.background= element_rect(color="white", size = 2, fill = "#f9f9f9"), 
        strip.text.y.left = element_text(angle = 0),
        strip.placement = "outside",
        strip.text.y = element_text(colour = "#091926", size = 12, face = "bold")) +
  xlab("Change: Pr(choice DMP)") +
  ggtitle("Tasks 2-5") +
  geom_point(size=2) +
  scale_color_manual(values = rep("#515151", 9))

combined_task <- ggarrange(c + theme(axis.ticks.y = element_blank(),
                                     plot.margin = margin(r = 1) ),
                           c1 + theme(axis.text.y = element_blank(),
                                      axis.ticks.y = element_blank(),
                                      axis.title.y = element_blank(),
                                      plot.margin = margin(l = 1)),
                           c2 + theme(axis.text.y = element_blank(),
                                      axis.ticks.y = element_blank(),
                                      axis.title.y = element_blank(),
                                      plot.margin = margin(l = 1)),
                           nrow = 1)

# Table E. 10: Test for carry-over effects
C29a <- head(amce_by_tasks [c("level", "estimate", "std.error")], 95L)
knitr::kable(C29a, format = "markdown")
C29b <- head(choice_amce_tasks_1 [c("level", "estimate", "std.error")], 10L)
knitr::kable(C29b, format = "markdown")
C29c <- head(choice_amce_tasks [c("level", "estimate", "std.error")], 10L)
knitr::kable(C29c, format = "markdown")

# Figure E. 11: Test for randomization / balancing
p_trust_politics <- plot(amce(cj_trust_politics, trust_politics ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav, id = ~ID))

politics <- amce(cj_trust_politics, trust_politics ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav, id = ~ID)

p_s1 <- p_trust_politics+theme(panel.border = element_rect(linetype = "blank"),
                               panel.background = element_blank(), 
                               axis.text.y = element_text(colour = "#515151", size = 12),
                               axis.text.x = element_text(colour = "#091926", size = 12),
                               axis.title = element_text(colour = "#091926", size = 12),
                               panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
                               axis.ticks = element_blank(), 
                               legend.position = "none",
                               panel.spacing.y = unit(1, "lines"),
                               strip.background= element_rect(color="white", size = 2, fill = "#f9f9f9"), 
                               strip.text.y.left = element_text(angle = 0),
                               strip.placement = "outside",
                               strip.text.y = element_text(colour = "#091926", size = 12, face = "bold")) +
  ggtitle("Trust in political institutions")+
  scale_color_manual(values = rep("#515151", 9))

p_trust_politicians <- plot(amce(cj_trust_politicians, trust_politicians ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav, id = ~ID))

politicians <- amce(cj_trust_politicians, trust_politicians ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav, id = ~ID)

p_s2 <- p_trust_politicians+theme(panel.border = element_rect(linetype = "blank"),
                                  panel.background = element_blank(), 
                                  axis.text.y = element_text(colour = "#515151", size = 12),
                                  axis.text.x = element_text(colour = "#091926", size = 12),
                                  axis.title = element_text(colour = "#091926", size = 12),
                                  panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
                                  axis.ticks = element_blank(), 
                                  legend.position = "none",
                                  panel.spacing.y = unit(1, "lines"),
                                  strip.background= element_rect(color="white", size = 2, fill = "#f9f9f9"), 
                                  strip.text.y.left = element_text(angle = 0),
                                  strip.placement = "outside",
                                  strip.text.y = element_text(colour = "#091926", size = 12, face = "bold")) +
  ggtitle("Trust in politicians")+
  scale_color_manual(values = rep("#515151", 9))

p_trust_citizens <- plot(amce(cj_trust, trust_citizens ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav, id = ~ID))

citizens <- amce(cj_trust, trust_citizens ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav, id = ~ID)

p_s3 <- p_trust_citizens+theme(panel.border = element_rect(linetype = "blank"),
                               panel.background = element_blank(), 
                               axis.text.y = element_text(colour = "#515151", size = 12),
                               axis.text.x = element_text(colour = "#091926", size = 12),
                               axis.title = element_text(colour = "#091926", size = 12),
                               panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
                               axis.ticks = element_blank(), 
                               legend.position = "none",
                               panel.spacing.y = unit(1, "lines"),
                               strip.background= element_rect(color="white", size = 2, fill = "#f9f9f9"), 
                               strip.text.y.left = element_text(angle = 0),
                               strip.placement = "outside",
                               strip.text.y = element_text(colour = "#091926", size = 12, face = "bold")) +
  ggtitle("Trust in citizens")+
  scale_color_manual(values = rep("#515151", 9))

p_trust_dmp <- plot(amce(cj_trust, trust_dmp ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav, id = ~ID))

dmp <- amce(cj_trust, trust_dmp ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav, id = ~ID)

p_s4 <- p_trust_dmp+theme(panel.border = element_rect(linetype = "blank"),
                          panel.background = element_blank(), 
                          axis.text.y = element_text(colour = "#515151", size = 12),
                          axis.text.x = element_text(colour = "#091926", size = 12),
                          axis.title = element_text(colour = "#091926", size = 12),
                          panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
                          axis.ticks = element_blank(), 
                          legend.position = "none",
                          panel.spacing.y = unit(1, "lines"),
                          strip.background= element_rect(color="white", size = 2, fill = "#f9f9f9"), 
                          strip.text.y.left = element_text(angle = 0),
                          strip.placement = "outside",
                          strip.text.y = element_text(colour = "#091926", size = 12, face = "bold")) +
  ggtitle("Trust in minipublics")+
  scale_color_manual(values = rep("#515151", 9))

balance <- ggarrange(p_s1 + theme(axis.ticks.y = element_blank(),
                                  plot.margin = margin(r = 1) ),
                     p_s2 + theme(axis.text.y = element_blank(),
                                  axis.ticks.y = element_blank(),
                                  axis.title.y = element_blank(),
                                  plot.margin = margin(l = 1)),
                     p_s4 + theme(axis.ticks.y = element_blank(),
                                  plot.margin = margin(l = 1)),
                     p_s3 + theme(axis.text.y = element_blank(),
                                  axis.ticks.y = element_blank(),
                                  axis.title.y = element_blank(),
                                  plot.margin = margin(l = 1)),
                     nrow = 2)

# Table E. 11: Test for randomization / balancing
C30a <- head(politics [c("level", "estimate", "std.error")], 95L)
knitr::kable(C30a, format = "markdown")
C30b <- head(politicians [c("level", "estimate", "std.error")], 19L)
knitr::kable(C30b, format = "markdown")
C30c <- head(politicians [c("level", "estimate", "std.error")], 19L)
knitr::kable(C30c, format = "markdown")
C30d <- head(dmp [c("level", "estimate", "std.error")], 19L)
knitr::kable(C30d, format = "markdown")

# Figure E. 12: Test for profile order effects
class(cj_all$profile)
cj_all$profile <- factor(cj_all$profile, labels = c("profile 1", "profile 2"))
class(cj_all$profile)

summary(cj_all$profile)

amce_profiles <- cj(cj_all, Mchoice, id = ~ID, estimate="amce_diff", feature_labels = feature_names, feature_order = c("outcome.fav", "authorization", "aim", "majority", "composition", "size", "recruitment", "initiative", "issue"), by = ~profile)

plot(amce_profiles, feature_headers = FALSE) +
  ggplot2::facet_wrap(~BY, ncol = 2) +
  theme(strip.background.x = element_blank(),
        strip.text.x = element_text(colour = "#091926", size = 14),
        strip.text.y = element_blank(),
        panel.spacing.y = unit(1, "lines"),
        panel.border = element_rect(linetype = "blank"),
        panel.background = element_blank(), 
        axis.text.y = element_text(colour = "#515151", size = 12),
        axis.text.x = element_text(colour = "#091926", size = 12),
        axis.title = element_text(colour = "#091926", size = 12),
        panel.grid.major = element_line(colour = "grey", linetype = "dotted"),
        axis.ticks = element_blank(), 
        legend.position = "none")+
  xlab("Change: Pr(choice DMP)") +
  geom_point(size=2)+
  scale_color_manual(values = rep("#515151",4)) 

# Table E. 12: Test for profile order effects
C31 <- head(amce_profiles [c("level", "estimate", "std.error")], 10L)
knitr::kable(C31, format = "markdown")

#### Multiple hypothesis testing ####

# SET WORKING DIRECTORY TO FOLDER CONTAINING REPLICATION FILES
# create working directory
working_directory <- ""
setwd(working_directory)

## Reproduce: 
##      - Figure E.13 - Table E.24
##      - Table E.13 - Table R.24  
##
## Output: 
##      Figure E.13: Estimates with corrections – authorization
##	    Table E.13: Estimates with corrections – authorization
##	    Figure E.14: Estimates with corrections – autonomy
##      Table E.14: Estimates with corrections – autonomy
##      Figure E.15: Estimates with corrections – additional provisions
##      Table E.15: Estimates with corrections – additional provisions
##      Figure E.16: Estimates with corrections – rest
##      Table E.16: Estimates with corrections – rest
##      Figure E.17: Estimates with corrections – People who have participated
##      Table E.17: Estimates with corrections – People who have participated
##      Figure E.18: Estimates with corrections – People who have not participated
##      Table E.18: Estimates with corrections – People who have not participated
##      Figure E.19: Estimates with corrections – People with high trust in minipublics
##      Table E.19: Estimates with corrections – People with high trust in minipublics
##      Figure E.20: Estimates with corrections – People with low trust in minipublics
##      Table E.20: Estimates with corrections – People with low trust in minipublics
##      Figure E.21: Estimates with corrections – People with high trust in citizens
##      Table E.21: Estimates with corrections – People with high trust in citizens
##      Figure E.22: Estimates with corrections – People with low trust in citizens
##      Table E.22: Estimates with corrections – People with low trust in citizens
##      Figure E.23: Estimates with corrections – People with high trust in citizens/minipublics but low trust in politics
##      Table E.23: Estimates with corrections – People with high trust in citizens/minipublics but low trust in politics
##      Figure E.24: Estimates with corrections – People who do not belong to the category high trust in citizens/minipublics but  low trust in politics
##      Table E.24: Estimates with corrections – People who do not belong to the category high trust in citizens/minipublics but  low trust in politics
##
## packages 
library(remotes)
library(cjoint)
library(ashr)
library(RColorBrewer) 
library(ggplot2)
library(plyr) 
library(haven) 
library(miceadds) 
library(DescTools) 
library(dplyr) 
library(varhandle) 
library(stringr) 
library(cregg)

# load in date
cj_all <- read_sav ("conjointdata_long_combined.sav") # choose correct path

## prepare conjoint attribute variables:
# as factors +
# labeling +
# rename attributes
cj_all$issue<-factor(cj_all$issue, labels = c("Emissions", "Immigration"))
cj_all$initiative <-factor(cj_all$initiative, labels = c("Think Tank/NGO", "Government"))
cj_all$recruitment <-factor(cj_all$recruitment, labels = c("Random selection", "Self-selection"))
cj_all$size <-factor(cj_all$size, labels = c("Small", "Large"))
cj_all$composition <-factor(cj_all$composition, labels = c("Citizens alone", "Mixed groups"))
cj_all$majority <-factor(cj_all$majority, labels = c("Narrow majority", "Clear majority"))
cj_all$aim <-factor(cj_all$aim, labels = c("Efficiency", "Consideration of all interests"))
cj_all$outcome.fav <- factor(cj_all$outcome.fav, labels = c("Preference mismatch","Preference match"))
cj_all$authorization <-factor(cj_all$authorization, labels = c("Binding decision",  "Recommendation to public officials", "Recommendation to a public referendum"))
cj_all$authorization <- relevel(cj_all$authorization, ref = "Recommendation to public officials")
cj_all$ID <- as.factor(cj_all$respondent_id)

feature_names <- c(
  'issue' = "Policy Issue",
  'initiative' = "Initiative",
  'recruitment' = "Recruitment",
  'size' = "Group Size",
  'composition' = "Group Composition",
  'majority' = "Degree of Consensus",
  'aim' = "Aim of the forum",
  'authorization' = "Authorization",
  'outcome.fav' = "Outcome Favorability"
)

## specifying conjoint models +
mod <- cjoint::amce(chosen ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav,
                    data=cj_all, cluster=TRUE, respondent.id="ID")

mod.sum <- summary.amce(mod)
amce.plot <- mod.sum$amce[2:4]

amce.plot <- amce.plot[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8),] # order it in the same way as the original paper


## get point estimates and se 
beta <- as.numeric(amce.plot$Estimate)
beta.se <- as.numeric(amce.plot$`Std. Err`)
names(beta) <- StrCap(amce.plot$Level)
beta.ref <- c("outcome.fav\n(ref: Preference mismatch)",
              "authorization\n (ref: Recommendation to public officials)",
              "aim\n (ref: Efficiency)",
              "majority\n (ref: Narrow majority)",
              "composition\n (ref: Citizens alone)",
              "size\n (ref: Small)",
              "recruitment\n (ref: Random selection)",
              "initiative\n (ref: Think Tank/NGO)",
              "issue\n (ref: Emissions)")

names(beta.se) <- amce.plot$Level

## get BH results
pVal <- mod.sum$amce[,6]
pVal <- pVal[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8)]
names(pVal) <- amce.plot$Level

p.bh <- p.adjust(pVal, "fdr")
p.bh <- ifelse(p.bh<=0.05, 1, 0)

## compute ASH results 

ashR.unif <- ash(beta, beta.se, mixcompdist = "uniform")
postBeta.unif <- ashR.unif$result[,9]
postBeta.sd.unif <- ashR.unif$result[,10]

ashR.norm <- ash(beta, beta.se, mixcompdist = "normal")
postBeta.norm <- ashR.norm$result[,9]
postBeta.sd.norm <- ashR.norm$result[,10]

## set threshold
alpha <- 0.05

test.num <- nrow(mod.sum$amce[2])

# Figure E.13: Estimates with corrections – authorization
# Table E.13: Estimates with corrections – authorization
beta_authorization <- beta[c(8:9)]
beta.se_authorization <- beta.se[c(8:9)]

postBeta.unif_authorization <- ashR.unif$result[c(8:9), 9]
postBeta.sd.unif_authorization <- ashR.unif$result[c(8:9), 10]

postBeta.norm_authorization <- ashR.norm$result[c(8:9), 9]
postBeta.sd.norm_authorization <- ashR.norm$result[c(8:9), 10]

p.bh_authorization <- p.bh[c(8:9)]

y.lab <- c("Binding decision", "Rec. to referendum")
y.length <- length(beta_authorization)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))
plot(NA, xlim = c(-0.08, 0.08), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")


abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")
segments(x0 = (beta_authorization - qnorm(1-alpha/2)*beta.se_authorization), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_authorization + qnorm(1-alpha/2)*beta.se_authorization), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_authorization - qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_authorization + qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_authorization - qnorm(0.975)*postBeta.sd.unif_authorization), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_authorization + qnorm(0.975)*postBeta.sd.unif_authorization), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_authorization - qnorm(0.975)*postBeta.sd.norm_authorization), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_authorization + qnorm(0.975)*postBeta.sd.norm_authorization), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_authorization, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_authorization, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_authorization, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_authorization, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_authorization+0.025)[p.bh_authorization==1], y.ind[y.ind.seq*p.bh_authorization]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = rep("Recommendation to \npublic officials", 2), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

# Figure E.14: Estimates with corrections – autonomy
# Table E.14: Estimates with corrections – autonomy

beta_autonomy <- beta[c(2, 5)]
beta.se_autonomy <- beta.se[c(2, 5)]

p.bh_autonomy <- p.bh[c(2, 5)]

postBeta.unif_autonomy <- ashR.unif$result[c(2, 5), 9]
postBeta.sd.unif_autonomy <- ashR.unif$result[c(2, 5), 10]

postBeta.norm_autonomy <- ashR.norm$result[c(2, 5), 9]
postBeta.sd.norm_autonomy <- ashR.norm$result[c(2, 5), 10]

y.lab <- c("Government", "Mixed groups")
y.length <- length(beta_autonomy)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.08, 0.08), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_autonomy - qnorm(1-alpha/2)*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_autonomy + qnorm(1-alpha/2)*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_autonomy - qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_autonomy + qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_autonomy - qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_autonomy + qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_autonomy - qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_autonomy + qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_autonomy, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_autonomy, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_autonomy, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_autonomy, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_autonomy+0.015)[p.bh_autonomy==1], y.ind[y.ind.seq*p.bh_autonomy]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = c("Think Tank / NGO", "Citizens alone"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

# Figure E.15: Estimates with corrections – additional provisions
# Table E.15: Estimates with corrections – additional provisions
beta_extra <- beta[c(3, 4, 6)]
beta.se_extra <- beta.se[c(3, 4, 6)]

p.bh_extra <- p.bh[c(3, 4, 6)]

postBeta.unif_extra <- ashR.unif$result[c(3, 4, 6), 9]
postBeta.sd.unif_extra <- ashR.unif$result[c(3, 4, 6), 10]

postBeta.norm_extra <- ashR.norm$result[c(3, 4, 6), 9]
postBeta.sd.norm_extra <- ashR.norm$result[c(3, 4, 6), 10]

y.lab <- c("Self-selection", "Large group", "Clear majority")
y.length <- length(beta_extra)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4, 5.5, length.out = 4),
           seq(7, 8.5, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(5,5,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.08, 0.08), 
     ylim = c(0.4, 8.5), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_extra - qnorm(1-alpha/2)*beta.se_extra), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_extra + qnorm(1-alpha/2)*beta.se_extra), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_extra - qnorm(1-alpha/(2*test.num))*beta.se_extra), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_extra + qnorm(1-alpha/(2*test.num))*beta.se_extra), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_extra - qnorm(0.975)*postBeta.sd.unif_extra), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_extra + qnorm(0.975)*postBeta.sd.unif_extra), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_extra - qnorm(0.975)*postBeta.sd.norm_extra), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_extra + qnorm(0.975)*postBeta.sd.norm_extra), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_extra, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_extra, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_extra, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_extra, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_extra+0.015)[p.bh_extra==1], y.ind[y.ind.seq*p.bh_extra]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.3, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.3, labels = c("Random selection", "Small groups", "Narrow majority"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 3), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 3), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 3), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 3), 
     las = 2, tick = F, cex.axis = 0.6)

# Figure E.16: Estimates with corrections – rest
# Table E.16: Estimates with corrections – rest
beta_rest <- beta[c(7, 1, 10)]
beta.se_rest <- beta.se[c(7, 1, 10)]

p.bh_rest <- p.bh[c(7, 1, 10)]

postBeta.unif_rest <- ashR.unif$result[c(7, 1, 10), 9]
postBeta.sd.unif_rest <- ashR.unif$result[c(7, 1, 10), 10]

postBeta.norm_rest <- ashR.norm$result[c(7, 1, 10), 9]
postBeta.sd.norm_rest <- ashR.norm$result[c(7, 1, 10), 10]

y.lab <- c("All interests", "Immigration", "Preference match")
y.length <- length(beta_rest)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4, 5.5, length.out = 4),
           seq(7, 8.5, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(5,5,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.08, 0.08), 
     ylim = c(0.4, 8.5), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_rest - qnorm(1-alpha/2)*beta.se_rest), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_rest + qnorm(1-alpha/2)*beta.se_rest), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_rest - qnorm(1-alpha/(2*test.num))*beta.se_rest), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_rest + qnorm(1-alpha/(2*test.num))*beta.se_rest), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_rest - qnorm(0.975)*postBeta.sd.unif_rest), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_rest + qnorm(0.975)*postBeta.sd.unif_rest), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_rest - qnorm(0.975)*postBeta.sd.norm_rest), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_rest + qnorm(0.975)*postBeta.sd.norm_rest), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_rest, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_rest, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_rest, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_rest, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_rest+0.015)[p.bh_rest==1], y.ind[y.ind.seq*p.bh_rest]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.3, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.3, labels = c("Efficiency", "Emissions", "Preference \nmissmatch"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 3), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 3), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 3), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 3), 
     las = 2, tick = F, cex.axis = 0.6)

## Participation
cj_participated <- cj_experience %>% filter(participated == 1)

## specifying conjoint models +
mod <- cjoint::amce(chosen ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav,
                    data=cj_participated, cluster=TRUE, respondent.id="ID")

mod.sum <- summary.amce(mod)
amce.plot <- mod.sum$amce[2:4]

amce.plot <- amce.plot[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8),] # order it in the same way as the original paper

## get point estimates and se 
beta <- as.numeric(amce.plot$Estimate)
beta.se <- as.numeric(amce.plot$`Std. Err`)
names(beta) <- StrCap(amce.plot$Level)
beta.ref <- c("outcome.fav\n(ref: Preference mismatch)",
              "authorization\n (ref: Recommendation to public officials)",
              "aim\n (ref: Efficiency)",
              "majority\n (ref: Narrow majority)",
              "composition\n (ref: Citizens alone)",
              "size\n (ref: Small)",
              "recruitment\n (ref: Random selection)",
              "initiative\n (ref: Think Tank/NGO)",
              "issue\n (ref: Emissions)")

names(beta.se) <- amce.plot$Level

## get BH results
pVal <- mod.sum$amce[,6]
pVal <- pVal[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8)]
names(pVal) <- amce.plot$Level

p.bh <- p.adjust(pVal, "fdr")
p.bh <- ifelse(p.bh<=0.05, 1, 0)

## compute ASH results 

ashR.unif <- ash(beta, beta.se, mixcompdist = "uniform")
postBeta.unif <- ashR.unif$result[,9]
postBeta.sd.unif <- ashR.unif$result[,10]

ashR.norm <- ash(beta, beta.se, mixcompdist = "normal")
postBeta.norm <- ashR.norm$result[,9]
postBeta.sd.norm <- ashR.norm$result[,10]


## set threshold
alpha <- 0.05

test.num <- nrow(mod.sum$amce[2])

# Figure E.17: Estimates with corrections – People who have participated
# Table E.17: Estimates with corrections – People who have participated
# authorization
beta_authorization <- beta[c(8:9)]
beta.se_authorization <- beta.se[c(8:9)]

postBeta.unif_authorization <- ashR.unif$result[c(8:9), 9]
postBeta.sd.unif_authorization <- ashR.unif$result[c(8:9), 10]

postBeta.norm_authorization <- ashR.norm$result[c(8:9), 9]
postBeta.sd.norm_authorization <- ashR.norm$result[c(8:9), 10]

p.bh_authorization <- p.bh[c(8:9)]

y.lab <- c("Binding decision", "Rec. to referendum")
y.length <- length(beta_authorization)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))
plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")


abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")
segments(x0 = (beta_authorization - qnorm(1-alpha/2)*beta.se_authorization), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_authorization + qnorm(1-alpha/2)*beta.se_authorization), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_authorization - qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_authorization + qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_authorization - qnorm(0.975)*postBeta.sd.unif_authorization), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_authorization + qnorm(0.975)*postBeta.sd.unif_authorization), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_authorization - qnorm(0.975)*postBeta.sd.norm_authorization), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_authorization + qnorm(0.975)*postBeta.sd.norm_authorization), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_authorization, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_authorization, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_authorization, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_authorization, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
#text((beta_authorization+0.025)[p.bh_authorization==1], y.ind[y.ind.seq*p.bh_authorization]+0.15, 
#     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = rep("Recommendation to \npublic officials", 2), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

# autonomy
beta_autonomy <- beta[c(2, 5)]
beta.se_autonomy <- beta.se[c(2, 5)]

p.bh_autonomy <- p.bh[c(2, 5)]

postBeta.unif_autonomy <- ashR.unif$result[c(2, 5), 9]
postBeta.sd.unif_autonomy <- ashR.unif$result[c(2, 5), 10]

postBeta.norm_autonomy <- ashR.norm$result[c(2, 5), 9]
postBeta.sd.norm_autonomy <- ashR.norm$result[c(2, 5), 10]

y.lab <- c("Government", "Mixed groups")
y.length <- length(beta_autonomy)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_autonomy - qnorm(1-alpha/2)*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_autonomy + qnorm(1-alpha/2)*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_autonomy - qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_autonomy + qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_autonomy - qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_autonomy + qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_autonomy - qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_autonomy + qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_autonomy, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_autonomy, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_autonomy, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_autonomy, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
#text((beta_autonomy+0.015)[p.bh_autonomy==1], y.ind[y.ind.seq*p.bh_autonomy]+0.15, 
#     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = c("Think Tank / NGO", "Citizens alone"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

## No participation
cj_notparticipated <- cj_experience %>% filter(participated == 0)

## specifying conjoint models +
mod <- cjoint::amce(chosen ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav,
                    data=cj_notparticipated, cluster=TRUE, respondent.id="ID")

mod.sum <- summary.amce(mod)
amce.plot <- mod.sum$amce[2:4]

amce.plot <- amce.plot[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8),] # order it in the same way as the original paper


## get point estimates and se 
beta <- as.numeric(amce.plot$Estimate)
beta.se <- as.numeric(amce.plot$`Std. Err`)
names(beta) <- StrCap(amce.plot$Level)
beta.ref <- c("outcome.fav\n(ref: Preference mismatch)",
              "authorization\n (ref: Recommendation to public officials)",
              "aim\n (ref: Efficiency)",
              "majority\n (ref: Narrow majority)",
              "composition\n (ref: Citizens alone)",
              "size\n (ref: Small)",
              "recruitment\n (ref: Random selection)",
              "initiative\n (ref: Think Tank/NGO)",
              "issue\n (ref: Emissions)")

names(beta.se) <- amce.plot$Level


## get BH results
pVal <- mod.sum$amce[,6]
pVal <- pVal[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8)]
names(pVal) <- amce.plot$Level

p.bh <- p.adjust(pVal, "fdr")
p.bh <- ifelse(p.bh<=0.05, 1, 0)

## compute ASH results 

ashR.unif <- ash(beta, beta.se, mixcompdist = "uniform")
postBeta.unif <- ashR.unif$result[,9]
postBeta.sd.unif <- ashR.unif$result[,10]

ashR.norm <- ash(beta, beta.se, mixcompdist = "normal")
postBeta.norm <- ashR.norm$result[,9]
postBeta.sd.norm <- ashR.norm$result[,10]


## set threshold
alpha <- 0.05

test.num <- nrow(mod.sum$amce[2])

# Figure E.18: Estimates with corrections – People who have not participated
# Table E.18: Estimates with corrections – People who have not participated
# authorization
beta_authorization <- beta[c(8:9)]
beta.se_authorization <- beta.se[c(8:9)]

postBeta.unif_authorization <- ashR.unif$result[c(8:9), 9]
postBeta.sd.unif_authorization <- ashR.unif$result[c(8:9), 10]

postBeta.norm_authorization <- ashR.norm$result[c(8:9), 9]
postBeta.sd.norm_authorization <- ashR.norm$result[c(8:9), 10]

p.bh_authorization <- p.bh[c(8:9)]

y.lab <- c("Binding decision", "Rec. to referendum")
y.length <- length(beta_authorization)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))
plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")


abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")
segments(x0 = (beta_authorization - qnorm(1-alpha/2)*beta.se_authorization), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_authorization + qnorm(1-alpha/2)*beta.se_authorization), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_authorization - qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_authorization + qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_authorization - qnorm(0.975)*postBeta.sd.unif_authorization), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_authorization + qnorm(0.975)*postBeta.sd.unif_authorization), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_authorization - qnorm(0.975)*postBeta.sd.norm_authorization), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_authorization + qnorm(0.975)*postBeta.sd.norm_authorization), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_authorization, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_authorization, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_authorization, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_authorization, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_authorization+0.025)[p.bh_authorization==1], y.ind[y.ind.seq*p.bh_authorization]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = rep("Recommendation to \npublic officials", 2), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

## autonomy
beta_autonomy <- beta[c(2, 5)]
beta.se_autonomy <- beta.se[c(2, 5)]

p.bh_autonomy <- p.bh[c(2, 5)]

postBeta.unif_autonomy <- ashR.unif$result[c(2, 5), 9]
postBeta.sd.unif_autonomy <- ashR.unif$result[c(2, 5), 10]

postBeta.norm_autonomy <- ashR.norm$result[c(2, 5), 9]
postBeta.sd.norm_autonomy <- ashR.norm$result[c(2, 5), 10]

y.lab <- c("Government", "Mixed groups")
y.length <- length(beta_autonomy)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_autonomy - qnorm(1-alpha/2)*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_autonomy + qnorm(1-alpha/2)*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_autonomy - qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_autonomy + qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_autonomy - qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_autonomy + qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_autonomy - qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_autonomy + qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_autonomy, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_autonomy, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_autonomy, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_autonomy, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_autonomy+0.015)[p.bh_autonomy==1], y.ind[y.ind.seq*p.bh_autonomy]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = c("Think Tank / NGO", "Citizens alone"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

# high trust in minipublics
summary(cj_trust$dmpsplit)
cj_trustdmp_high <- cj_trust %>% filter(dmpsplit == "High trust")

## specifying conjoint models
mod <- cjoint::amce(chosen ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav,
                    data=cj_trustdmp_high, cluster=TRUE, respondent.id="ID")

mod.sum <- summary.amce(mod)

amce.plot <- mod.sum$amce[2:4]

amce.plot <- amce.plot[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8),] # order it in the same way as the original paper


## get point estimates and se 
beta <- as.numeric(amce.plot$Estimate)
beta.se <- as.numeric(amce.plot$`Std. Err`)
names(beta) <- StrCap(amce.plot$Level)
beta.ref <- c("outcome.fav\n(ref: Preference mismatch)",
              "authorization\n (ref: Recommendation to public officials)",
              "aim\n (ref: Efficiency)",
              "majority\n (ref: Narrow majority)",
              "composition\n (ref: Citizens alone)",
              "size\n (ref: Small)",
              "recruitment\n (ref: Random selection)",
              "initiative\n (ref: Think Tank/NGO)",
              "issue\n (ref: Emissions)")

names(beta.se) <- amce.plot$Level


## get BH results
pVal <- mod.sum$amce[,6]
pVal <- pVal[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8)]
names(pVal) <- amce.plot$Level

p.bh <- p.adjust(pVal, "fdr")
p.bh <- ifelse(p.bh<=0.05, 1, 0)

## compute ASH results 

ashR.unif <- ash(beta, beta.se, mixcompdist = "uniform")
postBeta.unif <- ashR.unif$result[,9]
postBeta.sd.unif <- ashR.unif$result[,10]

ashR.norm <- ash(beta, beta.se, mixcompdist = "normal")
postBeta.norm <- ashR.norm$result[,9]
postBeta.sd.norm <- ashR.norm$result[,10]


## set threshold
alpha <- 0.05

test.num <- nrow(mod.sum$amce[2])

# Figure E.19: Estimates with corrections – People with high trust in minipublics
# Table E.19: Estimates with corrections – People with high trust in minipublics
## authorization
beta_authorization <- beta[c(8:9)]
beta.se_authorization <- beta.se[c(8:9)]

postBeta.unif_authorization <- ashR.unif$result[c(8:9), 9]
postBeta.sd.unif_authorization <- ashR.unif$result[c(8:9), 10]

postBeta.norm_authorization <- ashR.norm$result[c(8:9), 9]
postBeta.sd.norm_authorization <- ashR.norm$result[c(8:9), 10]

p.bh_authorization <- p.bh[c(8:9)]

y.lab <- c("Binding decision", "Rec. to referendum")
y.length <- length(beta_authorization)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))
plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")


abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")
segments(x0 = (beta_authorization - qnorm(1-alpha/2)*beta.se_authorization), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_authorization + qnorm(1-alpha/2)*beta.se_authorization), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_authorization - qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_authorization + qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_authorization - qnorm(0.975)*postBeta.sd.unif_authorization), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_authorization + qnorm(0.975)*postBeta.sd.unif_authorization), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_authorization - qnorm(0.975)*postBeta.sd.norm_authorization), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_authorization + qnorm(0.975)*postBeta.sd.norm_authorization), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_authorization, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_authorization, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_authorization, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_authorization, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_authorization+0.025)[p.bh_authorization==1], y.ind[y.ind.seq*p.bh_authorization]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = rep("Recommendation to \npublic officials", 2), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

## autonomy
beta_autonomy <- beta[c(2, 5)]
beta.se_autonomy <- beta.se[c(2, 5)]

p.bh_autonomy <- p.bh[c(2, 5)]

postBeta.unif_autonomy <- ashR.unif$result[c(2, 5), 9]
postBeta.sd.unif_autonomy <- ashR.unif$result[c(2, 5), 10]

postBeta.norm_autonomy <- ashR.norm$result[c(2, 5), 9]
postBeta.sd.norm_autonomy <- ashR.norm$result[c(2, 5), 10]

y.lab <- c("Government", "Mixed groups")
y.length <- length(beta_autonomy)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_autonomy - qnorm(1-alpha/2)*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_autonomy + qnorm(1-alpha/2)*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_autonomy - qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_autonomy + qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_autonomy - qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_autonomy + qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_autonomy - qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_autonomy + qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_autonomy, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_autonomy, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_autonomy, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_autonomy, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_autonomy+0.015)[p.bh_autonomy==1], y.ind[y.ind.seq*p.bh_autonomy]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = c("Think Tank / NGO", "Citizens alone"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)


# low trust in minipublics
summary(cj_trust$dmpsplit)
cj_trustdmp_low <- cj_trust %>% filter(dmpsplit == "Low trust")

## specifying conjoint models
mod <- cjoint::amce(chosen ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav,
                    data=cj_trustdmp_low, cluster=TRUE, respondent.id="ID")

mod.sum <- summary.amce(mod)

amce.plot <- mod.sum$amce[2:4]

amce.plot <- amce.plot[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8),] # order it in the same way as the original paper


## get point estimates and se 
beta <- as.numeric(amce.plot$Estimate)
beta.se <- as.numeric(amce.plot$`Std. Err`)
names(beta) <- StrCap(amce.plot$Level)
beta.ref <- c("outcome.fav\n(ref: Preference mismatch)",
              "authorization\n (ref: Recommendation to public officials)",
              "aim\n (ref: Efficiency)",
              "majority\n (ref: Narrow majority)",
              "composition\n (ref: Citizens alone)",
              "size\n (ref: Small)",
              "recruitment\n (ref: Random selection)",
              "initiative\n (ref: Think Tank/NGO)",
              "issue\n (ref: Emissions)")

names(beta.se) <- amce.plot$Level


## get BH results
pVal <- mod.sum$amce[,6]
pVal <- pVal[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8)]
names(pVal) <- amce.plot$Level

p.bh <- p.adjust(pVal, "fdr")
p.bh <- ifelse(p.bh<=0.05, 1, 0)

## compute ASH results 

ashR.unif <- ash(beta, beta.se, mixcompdist = "uniform")
postBeta.unif <- ashR.unif$result[,9]
postBeta.sd.unif <- ashR.unif$result[,10]

ashR.norm <- ash(beta, beta.se, mixcompdist = "normal")
postBeta.norm <- ashR.norm$result[,9]
postBeta.sd.norm <- ashR.norm$result[,10]


## set threshold
alpha <- 0.05

test.num <- nrow(mod.sum$amce[2])

# Figure E.20: Estimates with corrections – People with low trust in minipublics
# Table E.20: Estimates with corrections – People with low trust in minipublics
## authorization
beta_authorization <- beta[c(8:9)]
beta.se_authorization <- beta.se[c(8:9)]

postBeta.unif_authorization <- ashR.unif$result[c(8:9), 9]
postBeta.sd.unif_authorization <- ashR.unif$result[c(8:9), 10]

postBeta.norm_authorization <- ashR.norm$result[c(8:9), 9]
postBeta.sd.norm_authorization <- ashR.norm$result[c(8:9), 10]

p.bh_authorization <- p.bh[c(8:9)]

y.lab <- c("Binding decision", "Rec. to referendum")
y.length <- length(beta_authorization)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))
plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")


abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")
segments(x0 = (beta_authorization - qnorm(1-alpha/2)*beta.se_authorization), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_authorization + qnorm(1-alpha/2)*beta.se_authorization), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_authorization - qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_authorization + qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_authorization - qnorm(0.975)*postBeta.sd.unif_authorization), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_authorization + qnorm(0.975)*postBeta.sd.unif_authorization), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_authorization - qnorm(0.975)*postBeta.sd.norm_authorization), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_authorization + qnorm(0.975)*postBeta.sd.norm_authorization), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_authorization, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_authorization, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_authorization, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_authorization, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_authorization+0.025)[p.bh_authorization==1], y.ind[y.ind.seq*p.bh_authorization]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = rep("Recommendation to \npublic officials", 2), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

## autonomy
beta_autonomy <- beta[c(2, 5)]
beta.se_autonomy <- beta.se[c(2, 5)]

p.bh_autonomy <- p.bh[c(2, 5)]

postBeta.unif_autonomy <- ashR.unif$result[c(2, 5), 9]
postBeta.sd.unif_autonomy <- ashR.unif$result[c(2, 5), 10]

postBeta.norm_autonomy <- ashR.norm$result[c(2, 5), 9]
postBeta.sd.norm_autonomy <- ashR.norm$result[c(2, 5), 10]

y.lab <- c("Government", "Mixed groups")
y.length <- length(beta_autonomy)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_autonomy - qnorm(1-alpha/2)*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_autonomy + qnorm(1-alpha/2)*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_autonomy - qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_autonomy + qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_autonomy - qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_autonomy + qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_autonomy - qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_autonomy + qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_autonomy, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_autonomy, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_autonomy, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_autonomy, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_autonomy+0.015)[p.bh_autonomy==1], y.ind[y.ind.seq*p.bh_autonomy]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = c("Think Tank / NGO", "Citizens alone"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)


## trust in other citizens
cj_trustcit_high <- cj_trust %>% filter(citizensplit == "High trust")

## specifying conjoint models +
mod <- cjoint::amce(chosen ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav,
                    data=cj_trustcit_high, cluster=TRUE, respondent.id="ID")

mod.sum <- summary.amce(mod)
amce.plot <- mod.sum$amce[2:4]

amce.plot <- amce.plot[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8),] # order it in the same way as the original paper


## get point estimates and se 
beta <- as.numeric(amce.plot$Estimate)
beta.se <- as.numeric(amce.plot$`Std. Err`)
names(beta) <- StrCap(amce.plot$Level)
beta.ref <- c("outcome.fav\n(ref: Preference mismatch)",
              "authorization\n (ref: Recommendation to public officials)",
              "aim\n (ref: Efficiency)",
              "majority\n (ref: Narrow majority)",
              "composition\n (ref: Citizens alone)",
              "size\n (ref: Small)",
              "recruitment\n (ref: Random selection)",
              "initiative\n (ref: Think Tank/NGO)",
              "issue\n (ref: Emissions)")

names(beta.se) <- amce.plot$Level


## get BH results
pVal <- mod.sum$amce[,6]
pVal <- pVal[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8)]
names(pVal) <- amce.plot$Level

p.bh <- p.adjust(pVal, "fdr")
p.bh <- ifelse(p.bh<=0.05, 1, 0)

## compute ASH results 

ashR.unif <- ash(beta, beta.se, mixcompdist = "uniform")
postBeta.unif <- ashR.unif$result[,9]
postBeta.sd.unif <- ashR.unif$result[,10]

ashR.norm <- ash(beta, beta.se, mixcompdist = "normal")
postBeta.norm <- ashR.norm$result[,9]
postBeta.sd.norm <- ashR.norm$result[,10]


## set threshold
alpha <- 0.05

test.num <- nrow(mod.sum$amce[2])

# Figure E.21: Estimates with corrections – People with high trust in citizens
# Table E.21: Estimates with corrections – People with high trust in citizens
## authorization
beta_authorization <- beta[c(8:9)]
beta.se_authorization <- beta.se[c(8:9)]

postBeta.unif_authorization <- ashR.unif$result[c(8:9), 9]
postBeta.sd.unif_authorization <- ashR.unif$result[c(8:9), 10]

postBeta.norm_authorization <- ashR.norm$result[c(8:9), 9]
postBeta.sd.norm_authorization <- ashR.norm$result[c(8:9), 10]

p.bh_authorization <- p.bh[c(8:9)]

y.lab <- c("Binding decision", "Rec. to referendum")
y.length <- length(beta_authorization)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))
plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")


abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")
segments(x0 = (beta_authorization - qnorm(1-alpha/2)*beta.se_authorization), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_authorization + qnorm(1-alpha/2)*beta.se_authorization), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_authorization - qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_authorization + qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_authorization - qnorm(0.975)*postBeta.sd.unif_authorization), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_authorization + qnorm(0.975)*postBeta.sd.unif_authorization), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_authorization - qnorm(0.975)*postBeta.sd.norm_authorization), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_authorization + qnorm(0.975)*postBeta.sd.norm_authorization), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_authorization, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_authorization, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_authorization, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_authorization, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_authorization+0.025)[p.bh_authorization==1], y.ind[y.ind.seq*p.bh_authorization]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = rep("Recommendation to \npublic officials", 2), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

## autonomy
beta_autonomy <- beta[c(2, 5)]
beta.se_autonomy <- beta.se[c(2, 5)]

p.bh_autonomy <- p.bh[c(2, 5)]

postBeta.unif_autonomy <- ashR.unif$result[c(2, 5), 9]
postBeta.sd.unif_autonomy <- ashR.unif$result[c(2, 5), 10]

postBeta.norm_autonomy <- ashR.norm$result[c(2, 5), 9]
postBeta.sd.norm_autonomy <- ashR.norm$result[c(2, 5), 10]

y.lab <- c("Government", "Mixed groups")
y.length <- length(beta_autonomy)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_autonomy - qnorm(1-alpha/2)*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_autonomy + qnorm(1-alpha/2)*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_autonomy - qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_autonomy + qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_autonomy - qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_autonomy + qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_autonomy - qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_autonomy + qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_autonomy, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_autonomy, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_autonomy, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_autonomy, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_autonomy+0.015)[p.bh_autonomy==1], y.ind[y.ind.seq*p.bh_autonomy]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = c("Think Tank / NGO", "Citizens alone"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)


# low trust in other citizens
cj_trustcit_low <- cj_trust %>% filter(citizensplit == "Low trust")

## specifying conjoint models
mod <- cjoint::amce(chosen ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav,
                    data=cj_trustcit_low, cluster=TRUE, respondent.id="ID")

mod.sum <- summary.amce(mod)

amce.plot <- mod.sum$amce[2:4]

amce.plot <- amce.plot[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8),] # order it in the same way as the original paper


## get point estimates and se 
beta <- as.numeric(amce.plot$Estimate)
beta.se <- as.numeric(amce.plot$`Std. Err`)
names(beta) <- StrCap(amce.plot$Level)
beta.ref <- c("outcome.fav\n(ref: Preference mismatch)",
              "authorization\n (ref: Recommendation to public officials)",
              "aim\n (ref: Efficiency)",
              "majority\n (ref: Narrow majority)",
              "composition\n (ref: Citizens alone)",
              "size\n (ref: Small)",
              "recruitment\n (ref: Random selection)",
              "initiative\n (ref: Think Tank/NGO)",
              "issue\n (ref: Emissions)")

names(beta.se) <- amce.plot$Level


## get BH results
pVal <- mod.sum$amce[,6]
pVal <- pVal[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8)]
names(pVal) <- amce.plot$Level

p.bh <- p.adjust(pVal, "fdr")
p.bh <- ifelse(p.bh<=0.05, 1, 0)

## compute ASH results 

ashR.unif <- ash(beta, beta.se, mixcompdist = "uniform")
postBeta.unif <- ashR.unif$result[,9]
postBeta.sd.unif <- ashR.unif$result[,10]

ashR.norm <- ash(beta, beta.se, mixcompdist = "normal")
postBeta.norm <- ashR.norm$result[,9]
postBeta.sd.norm <- ashR.norm$result[,10]


## set threshold
alpha <- 0.05

test.num <- nrow(mod.sum$amce[2])

# Figure E.22: Estimates with corrections – People with low trust in citizens
# Table E.22: Estimates with corrections – People with low trust in citizens
## authorization
beta_authorization <- beta[c(8:9)]
beta.se_authorization <- beta.se[c(8:9)]

postBeta.unif_authorization <- ashR.unif$result[c(8:9), 9]
postBeta.sd.unif_authorization <- ashR.unif$result[c(8:9), 10]

postBeta.norm_authorization <- ashR.norm$result[c(8:9), 9]
postBeta.sd.norm_authorization <- ashR.norm$result[c(8:9), 10]

p.bh_authorization <- p.bh[c(8:9)]

y.lab <- c("Binding decision", "Rec. to referendum")
y.length <- length(beta_authorization)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))
plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")


abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")
segments(x0 = (beta_authorization - qnorm(1-alpha/2)*beta.se_authorization), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_authorization + qnorm(1-alpha/2)*beta.se_authorization), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_authorization - qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_authorization + qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_authorization - qnorm(0.975)*postBeta.sd.unif_authorization), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_authorization + qnorm(0.975)*postBeta.sd.unif_authorization), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_authorization - qnorm(0.975)*postBeta.sd.norm_authorization), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_authorization + qnorm(0.975)*postBeta.sd.norm_authorization), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_authorization, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_authorization, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_authorization, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_authorization, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_authorization+0.025)[p.bh_authorization==1], y.ind[y.ind.seq*p.bh_authorization]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = rep("Recommendation to \npublic officials", 2), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

## autonomy
beta_autonomy <- beta[c(2, 5)]
beta.se_autonomy <- beta.se[c(2, 5)]

p.bh_autonomy <- p.bh[c(2, 5)]

postBeta.unif_autonomy <- ashR.unif$result[c(2, 5), 9]
postBeta.sd.unif_autonomy <- ashR.unif$result[c(2, 5), 10]

postBeta.norm_autonomy <- ashR.norm$result[c(2, 5), 9]
postBeta.sd.norm_autonomy <- ashR.norm$result[c(2, 5), 10]

y.lab <- c("Government", "Mixed groups")
y.length <- length(beta_autonomy)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_autonomy - qnorm(1-alpha/2)*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_autonomy + qnorm(1-alpha/2)*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_autonomy - qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_autonomy + qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_autonomy - qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_autonomy + qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_autonomy - qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_autonomy + qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_autonomy, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_autonomy, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_autonomy, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_autonomy, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_autonomy+0.015)[p.bh_autonomy==1], y.ind[y.ind.seq*p.bh_autonomy]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = c("Think Tank / NGO", "Citizens alone"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)


# low trust politics / high trust citizens
summary(cj_trust_NA$select_trust)
cj_trust_NA$select_trust <- as.factor(cj_trust_NA$select_trust)
cj_poltrust_low_soctrust_high<- cj_trust_NA %>% filter(select_trust == 1)

## specifying conjoint models
mod <- cjoint::amce(chosen ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav,
                    data=cj_poltrust_low_soctrust_high, cluster=TRUE, respondent.id="ID")

mod.sum <- summary.amce(mod)

amce.plot <- mod.sum$amce[2:4]

amce.plot <- amce.plot[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8),] # order it in the same way as the original paper


## get point estimates and se 
beta <- as.numeric(amce.plot$Estimate)
beta.se <- as.numeric(amce.plot$`Std. Err`)
names(beta) <- StrCap(amce.plot$Level)
beta.ref <- c("outcome.fav\n(ref: Preference mismatch)",
              "authorization\n (ref: Recommendation to public officials)",
              "aim\n (ref: Efficiency)",
              "majority\n (ref: Narrow majority)",
              "composition\n (ref: Citizens alone)",
              "size\n (ref: Small)",
              "recruitment\n (ref: Random selection)",
              "initiative\n (ref: Think Tank/NGO)",
              "issue\n (ref: Emissions)")

names(beta.se) <- amce.plot$Level


## get BH results
pVal <- mod.sum$amce[,6]
pVal <- pVal[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8)]
names(pVal) <- amce.plot$Level

p.bh <- p.adjust(pVal, "fdr")
p.bh <- ifelse(p.bh<=0.05, 1, 0)

## compute ASH results 

ashR.unif <- ash(beta, beta.se, mixcompdist = "uniform")
postBeta.unif <- ashR.unif$result[,9]
postBeta.sd.unif <- ashR.unif$result[,10]

ashR.norm <- ash(beta, beta.se, mixcompdist = "normal")
postBeta.norm <- ashR.norm$result[,9]
postBeta.sd.norm <- ashR.norm$result[,10]


## set threshold
alpha <- 0.05

test.num <- nrow(mod.sum$amce[2])

# Figure E.23: Estimates with corrections – People with high trust in citizens/minipublics but low trust in politics
# Table E.23: Estimates with corrections – People with high trust in citizens/minipublics but low trust in politics
## authorization
beta_authorization <- beta[c(8:9)]
beta.se_authorization <- beta.se[c(8:9)]

postBeta.unif_authorization <- ashR.unif$result[c(8:9), 9]
postBeta.sd.unif_authorization <- ashR.unif$result[c(8:9), 10]

postBeta.norm_authorization <- ashR.norm$result[c(8:9), 9]
postBeta.sd.norm_authorization <- ashR.norm$result[c(8:9), 10]

p.bh_authorization <- p.bh[c(8:9)]

y.lab <- c("Binding decision", "Rec. to referendum")
y.length <- length(beta_authorization)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))
plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")


abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")
segments(x0 = (beta_authorization - qnorm(1-alpha/2)*beta.se_authorization), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_authorization + qnorm(1-alpha/2)*beta.se_authorization), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_authorization - qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_authorization + qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_authorization - qnorm(0.975)*postBeta.sd.unif_authorization), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_authorization + qnorm(0.975)*postBeta.sd.unif_authorization), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_authorization - qnorm(0.975)*postBeta.sd.norm_authorization), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_authorization + qnorm(0.975)*postBeta.sd.norm_authorization), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_authorization, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_authorization, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_authorization, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_authorization, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_authorization+0.025)[p.bh_authorization==1], y.ind[y.ind.seq*p.bh_authorization]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = rep("Recommendation to \npublic officials", 2), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

## autonomy
beta_autonomy <- beta[c(2, 5)]
beta.se_autonomy <- beta.se[c(2, 5)]

p.bh_autonomy <- p.bh[c(2, 5)]

postBeta.unif_autonomy <- ashR.unif$result[c(2, 5), 9]
postBeta.sd.unif_autonomy <- ashR.unif$result[c(2, 5), 10]

postBeta.norm_autonomy <- ashR.norm$result[c(2, 5), 9]
postBeta.sd.norm_autonomy <- ashR.norm$result[c(2, 5), 10]

y.lab <- c("Government", "Mixed groups")
y.length <- length(beta_autonomy)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_autonomy - qnorm(1-alpha/2)*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_autonomy + qnorm(1-alpha/2)*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_autonomy - qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_autonomy + qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_autonomy - qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_autonomy + qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_autonomy - qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_autonomy + qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_autonomy, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_autonomy, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_autonomy, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_autonomy, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_autonomy+0.015)[p.bh_autonomy==1], y.ind[y.ind.seq*p.bh_autonomy]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = c("Think Tank / NGO", "Citizens alone"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)


# trust rest (high political / low social)
summary(cj_trust_NA$select_trust)
cj_trust_NA$select_trust <- as.factor(cj_trust_NA$select_trust)
cj_poltrust_soctrust_rest<- cj_trust_NA %>% filter(select_trust == 0)

## specifying conjoint models
mod <- cjoint::amce(chosen ~ issue + initiative + recruitment + size + composition + majority + aim + authorization + outcome.fav,
                    data=cj_poltrust_soctrust_rest, cluster=TRUE, respondent.id="ID")

mod.sum <- summary.amce(mod)

amce.plot <- mod.sum$amce[2:4]

amce.plot <- amce.plot[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8),] # order it in the same way as the original paper


## get point estimates and se 
beta <- as.numeric(amce.plot$Estimate)
beta.se <- as.numeric(amce.plot$`Std. Err`)
names(beta) <- StrCap(amce.plot$Level)
beta.ref <- c("outcome.fav\n(ref: Preference mismatch)",
              "authorization\n (ref: Recommendation to public officials)",
              "aim\n (ref: Efficiency)",
              "majority\n (ref: Narrow majority)",
              "composition\n (ref: Citizens alone)",
              "size\n (ref: Small)",
              "recruitment\n (ref: Random selection)",
              "initiative\n (ref: Think Tank/NGO)",
              "issue\n (ref: Emissions)")

names(beta.se) <- amce.plot$Level


## get BH results
pVal <- mod.sum$amce[,6]
pVal <- pVal[c(6, 5, 9, 10, 4, 7, 1, 2:3, 8)]
names(pVal) <- amce.plot$Level

p.bh <- p.adjust(pVal, "fdr")
p.bh <- ifelse(p.bh<=0.05, 1, 0)

## compute ASH results 

ashR.unif <- ash(beta, beta.se, mixcompdist = "uniform")
postBeta.unif <- ashR.unif$result[,9]
postBeta.sd.unif <- ashR.unif$result[,10]

ashR.norm <- ash(beta, beta.se, mixcompdist = "normal")
postBeta.norm <- ashR.norm$result[,9]
postBeta.sd.norm <- ashR.norm$result[,10]


## set threshold
alpha <- 0.05

test.num <- nrow(mod.sum$amce[2])

##      Figure E.24: Estimates with corrections – People who do not belong to the category high trust in citizens/minipublics but  low trust in politics
##      Table E.24: Estimates with corrections – People who do not belong to the category high trust in citizens/minipublics but  low trust in politics
## authorization
beta_authorization <- beta[c(8:9)]
beta.se_authorization <- beta.se[c(8:9)]

postBeta.unif_authorization <- ashR.unif$result[c(8:9), 9]
postBeta.sd.unif_authorization <- ashR.unif$result[c(8:9), 10]

postBeta.norm_authorization <- ashR.norm$result[c(8:9), 9]
postBeta.sd.norm_authorization <- ashR.norm$result[c(8:9), 10]

p.bh_authorization <- p.bh[c(8:9)]

y.lab <- c("Binding decision", "Rec. to referendum")
y.length <- length(beta_authorization)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))
plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")


abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")
segments(x0 = (beta_authorization - qnorm(1-alpha/2)*beta.se_authorization), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_authorization + qnorm(1-alpha/2)*beta.se_authorization), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_authorization - qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_authorization + qnorm(1-alpha/(2*test.num))*beta.se_authorization), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_authorization - qnorm(0.975)*postBeta.sd.unif_authorization), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_authorization + qnorm(0.975)*postBeta.sd.unif_authorization), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_authorization - qnorm(0.975)*postBeta.sd.norm_authorization), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_authorization + qnorm(0.975)*postBeta.sd.norm_authorization), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_authorization, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_authorization, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_authorization, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_authorization, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_authorization+0.025)[p.bh_authorization==1], y.ind[y.ind.seq*p.bh_authorization]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = rep("Recommendation to \npublic officials", 2), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)

## autonomy
beta_autonomy <- beta[c(2, 5)]
beta.se_autonomy <- beta.se[c(2, 5)]

p.bh_autonomy <- p.bh[c(2, 5)]

postBeta.unif_autonomy <- ashR.unif$result[c(2, 5), 9]
postBeta.sd.unif_autonomy <- ashR.unif$result[c(2, 5), 10]

postBeta.norm_autonomy <- ashR.norm$result[c(2, 5), 9]
postBeta.sd.norm_autonomy <- ashR.norm$result[c(2, 5), 10]

y.lab <- c("Government", "Mixed groups")
y.length <- length(beta_autonomy)
y.ind <- c(seq(1.5, 3, length.out = 4), seq(4.5, 6, length.out = 4))
y.ind <- sort(y.ind, decreasing = T)
y.ind.seq <- seq(1,length(y.ind), 4)

par(mar = c(4,8,4,8))
par(mai= c(1,1.5,1,1.5))

plot(NA, xlim = c(-0.2, 0.2), 
     ylim = c(0.4, 6), 
     type = "n",
     yaxt = "n", 
     xlab = "", 
     cex.lab = 1.2,
     cex.axis = 1.5, 
     ylab = "", axes = F)
abline(v = 0, col = "gray")
abline(h = y.ind[y.ind.seq], lwd = 1, col = "grey80")

segments(x0 = (beta_autonomy - qnorm(1-alpha/2)*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq], 
         x1 = (beta_autonomy + qnorm(1-alpha/2)*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq], col = "black", lwd = 2)
segments(x0 = (beta_autonomy - qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y0 = y.ind[y.ind.seq+1], 
         x1 = (beta_autonomy + qnorm(1-alpha/(2*test.num))*beta.se_autonomy), 
         y1 = y.ind[y.ind.seq+1], col = "black", lwd = 2)
segments(x0 = (postBeta.unif_autonomy - qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y0 = y.ind[y.ind.seq+2], 
         x1 = (postBeta.unif_autonomy + qnorm(0.975)*postBeta.sd.unif_autonomy), 
         y1 = y.ind[y.ind.seq+2], col = "black", lwd = 2)
segments(x0 = (postBeta.norm_autonomy - qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y0 = y.ind[y.ind.seq+3], 
         x1 = (postBeta.norm_autonomy + qnorm(0.975)*postBeta.sd.norm_autonomy), 
         y1 = y.ind[y.ind.seq+3], col = "black", lwd = 2)
points(beta_autonomy, y.ind[y.ind.seq], cex = 1.4, pch = 19, 
       col = "black")
points(beta_autonomy, y.ind[y.ind.seq+1], cex = 1.4, pch = 1, 
       col = "black")
points(postBeta.unif_autonomy, y.ind[y.ind.seq+2], cex = 1.4, pch = 6)
points(postBeta.norm_autonomy, y.ind[y.ind.seq+3], cex = 1.4, pch = 2)
text((beta_autonomy+0.015)[p.bh_autonomy==1], y.ind[y.ind.seq*p.bh_autonomy]+0.15, 
     labels = bquote("BH"~symbol("\326")),  cex = 0.75)
box()

axis(1, at = seq(-0.2,0.2, by=0.01), labels = seq(-0.2,0.2, by=0.01), cex.axis = 1)
axis(side = 4, at = y.ind[y.ind.seq]+0.22, labels = y.lab, las = 2, tick = F, cex.axis = 0.8, font = 0.8)
axis(side = 2, at = y.ind[y.ind.seq]+0.22, labels = c("Think Tank / NGO", "Citizens alone"), font = 0.8,
     las = 2, tick = F, cex.axis = 0.8)
axis(side = 4, at = y.ind[y.ind.seq], labels = rep("no Corr.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+1], labels = rep("Bonf.", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+2], labels = rep("ash.Unif", 2), 
     las = 2, tick = F, cex.axis = 0.6)
axis(side = 4, at = y.ind[y.ind.seq+3], labels = rep("ash.Norm", 2), 
     las = 2, tick = F, cex.axis = 0.6)
